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ABSTRACT 


Electron  transport  properties  of  terahertz  (THz)  longitudinal  optical  (LO)-phonon 
quantum  cascade  (QC)  structures  were  modeled,  in  order  to  investigate  high  gain 
quantum  cascade  laser  (QCL)  structures.  A  new  structure,  a  step  well  QC  structure,  was 
proposed.  Under  such  an  arrangement,  there  are  three  main  energy  levels  within  the  step 
well,  where  the  transition  from  the  upper  state  to  the  middle  state  is  at  the  THz  radiative 
spacing  and  the  transition  from  the  middle  state  to  the  lower  state  is  at  or  near  the  LO- 
phonon  energy  (~  36  meV  in  GaAs).  Because  of  the  inherent  difficulties  in  using  rate 
equations  for  this  type  of  transport  analysis,  a  Monte  Carlo  simulation  was  developed. 
Step  well  injectors  were  modeled  and  shown  to  be  capable  of  high  injection  efficiencies 
(~  90%),  higher  than  previously  obtained.  Comparisons  to  conventional  square  well  LO- 
phonon  structures  are  made,  including  a  Monte  Carlo  analysis  of  a  high  power  THz  QCL. 
Interface  roughness  scattering  was  shown  to  be  significant  only  for  roughness  greater 
than  approximately  one  monolayer.  It  was  found  that  step  well  structures  are  capable  of 
high  gains  and  injection  efficiencies,  with  comparable  characteristics  to  other  square  well 
designs,  but  do  have  increased  scattering  from  the  upper  state  to  the  lower  states. 
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/?=  complex  propagation  constant 

r  =  confinement  factor 

A,  =  wavelength 

£  =  “bed  spring”  spinor 

t=  lifetime 

y/=  “pitch  fork”  envelope  wavefunction 
co  =  angular  frequency 

Constants  taken  from  the  Particle  Physics  Booklet,  PDG  (July  2004). 
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EXECUTIVE  SUMMARY 


Electron  transport  properties  of  terahertz  (THz)  longitudinal  optical  (LO)-phonon 
quantum  cascade  (QC)  structures  were  modeled,  in  order  to  investigate  and  better 
understand  the  underlying  physics  of  high  gain  quantum  cascade  laser  (QCL)  structures. 
A  new  structure,  a  step  well  QC  structure,  was  proposed  after  rate  equation  analysis 
showed  favorable  scattering  rate  lifetimes  for  maintaining  a  population  inversion.  Under 
such  an  arrangement,  there  are  three  main  energy  levels  within  the  step  well,  where  the 
transition  from  the  upper  state  to  the  middle  state  is  at  the  THz  radiative  spacing  and  the 
transition  from  the  middle  state  to  the  lower  state  is  at  or  near  the  LO-phonon  energy 
(~  36  meV  in  GaAs)  for  fast  depopulation.  The  middle  state  (upper  phonon  or  lower 
lasing  state)  is  a  single  energy  state,  contrasting  to  previous  LO-phonon  based  QCL 
designs  that  have  doublet  states.  For  long  wavelength  lasers,  this  may  be  important 
because  it  eliminates  the  possibility  for  unwanted  THz  absorption  that  could  otherwise 
occur  between  those  doublet  states.  Since  the  radiative  and  LO-phonon  transitions  are 
intrawell  in  nature,  high  oscillator  strengths  and  sub-picosecond  middle  state  lifetimes  are 
possible,  which  can  lead  to  increased  gain  in  the  active  region  provided  the  upper  state 
lifetime  and  injection  efficiency  are  maintained.  Because  of  the  inherent  difficulties  in 
using  rate  equations  for  this  type  of  transport  analysis,  a  Monte  Carlo  simulation  was 
developed.  Step  well  injectors  were  modeled  and  shown  to  be  capable  of  high  injection 
efficiencies  (~  90%),  which  is  higher  than  previously  obtained  in  other  LO-phonon 
structures  with  similar  wavelengths  of  operation  and  oscillator  strengths.  It  is  likely  that 
step  injectors  could  be  useful  in  approaches  that  do  not  arrange  all  three  energy  levels 
within  the  same  well.  Comparisons  to  conventional  square  well  LO-phonon  structures  are 
made,  including  a  Monte  Carlo  analysis  of  a  high  power  THz  QCL.  Interface  roughness 
scattering  was  shown  to  be  significant  only  for  roughness  greater  than  approximately  one 
monolayer.  It  was  found  that  step  well  structures  are  capable  of  high  gains  and  oscillator 
strengths,  high  injection  efficiencies,  with  comparable  characteristics  to  other  square  well 
designs,  but  do  have  increased  scattering  from  the  upper  state  to  the  lower  states. 
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I.  INTRODUCTION 


Terahertz  (THz)  sources  capable  of  operating  in  roughly  the  1  to  5  THz  region 
(300  to  60  //in)  or  far-infrared  region  are  of  considerable  interest  and  are  currently  an 
active  field  of  study.  The  THz  frequency  range  has  remained  one  of  the  least  developed 
spectral  regions,  in  part  due  to  the  lack  of  compact  coherent  sources.  Though,  in  the  last 
seven  years  or  so,  there  has  been  a  significant  increase  in  efforts  to  produce  compact  THz 
sources.  One  approach  to  providing  a  compact  source  within  the  THz  region,  is  to  use  a 
quantum  cascade  laser  (QCL).1  QCLs  are  solid  state  unipolar  devices  that  differ  from 
conventional  semiconductor  lasers  in  that  transitions  occur  all  within  the  conduction 
band,  rather  than  between  the  bandgap  of  the  conduction  and  valence  bands.  This  class  of 
semiconductor  lasers  allows  for  customization  of  the  wavelength,  by  utilizing  a  multiple 
quantum  well  active  region.  A  key  feature  of  a  QCL  is  that  it  consists  of  N  repeated 
sections  that  form  a  cascade.  Thus,  one  injected  electron  can  emit  many  photons,  which 
allows  for  differential  quantum  efficiencies  greater  than  unity  and  hence  higher  power 
output. 

QCLs  have  been  successfully  designed  for  operating  in  the  infrared  (IR)  region 
(for  wavelengths  from  ~  3  to  24  //m,  with  room  temperature  operation,  and  CW  operation 
of  hundreds  of  mW  of  power).-'  More  recently,  QCLs  have  been  designed  in  the  THz 
region  (~  0.84  to  5  THz,  up  to  ~  185  K,  250  mW  pulsed,  and  130  mW  CW  operation).8'12 
Though,  for  THz  QCLs  these  figures  are  not  from  the  same  device.  Temperature  as  a 
function  of  frequency  for  the  best  devices  has  approximately  followed  the  trend  T  ~ 
E radiativj kB  =  hcolkB,  and  several  of  the  low  frequency  designs  operate  with  the  assistance 
of  a  magnetic  field.  Alternative  compact  solid  state  sources,  such  as  lead-salt 
semiconductor  diode  lasers,  are  limited  by  the  bandgap  to  ~  30  jui n  for  the  longest 
wavelengths.14  Other  solid  state  devices  such  as  transistors,  Gunn  oscillators,  and 
Schottky  diodes  multipliers,  have  at  best  achieved  low  juW  power  levels  in  the  THz 
region.15'18  Gas  lasers  generally  have  limited  lasing  frequencies  and  are  usually  fairly 
bulky.  For  these  reasons,  and  because  of  the  lack  of  intersubband  materials,  QCLs 
provide  a  good  means  for  solid  state  THz  sources. 
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Some  applications  for  the  terahertz  spectrum  include  astrophysical  science, 
medical  science,  THz  spectroscopy,  various  security  applications  such  as  THz  imaging 
systems,  chemical  gas  sensing,  and  agent  detection.15,19'31  Because  THz  radiation  is 
nonionizing,  it  is  attractive  for  security  imaging  applications.  Many  materials  have 
absorption  bands  in  the  THz  region  of  the  electromagnetic  spectrum,  therefore  this  may 
be  used  for  the  characterization  of  explosive  materials.  Imaging  in  this  region  can  allow 
for  characterization  of  materials  as  they  are  brought  through  a  security  screen,  as  well  as 
detect  and  identify  materials.  THz  waves  penetrate  dielectrics,  such  as  clothes  and  plastic, 
but  are  blocked  by  metallic  objects,  which  affect  their  screening  and  detection 
capabilities. 

A.  PREVIOUS  THz  QC  STRUCTURES 

THz  QCL  designs  have  used  two  main  quantum  cascade  (QC)  structure 
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architectures,  broadly  categorized  as  miniband  ’  '  ’  and  longitudinal  optical  (LO)- 
phonon,  ’  ’  using  coupled  square  quantum  wells  (square  when  unbiased).  Although 
sometimes  hybrid  structures  have  been  used  with  some  overlap  of  these  two  approaches, 
usually  designs  can  be  classified  as  belonging  to  one  of  these  two  classes. 

In  miniband  designs,  very  closely  spaced  energy  levels  are  formed  which  create 
minibands.  Radiative  transitions  take  place  between  these  minibands  or  between  a 
relatively  isolated  state  and  a  miniband.  The  first  THz  QCL  which  used  a  chirped 
superlattice,1  is  an  example  of  a  miniband  structure.  Other  miniband  designs  have 
included  a  bound  to  continuum  (Figure  I.A.l(a))  and  a  hybrid  bound  to  continuum  with 
optical-phonon  scattering  approaches.  "  Both  featured  minibands  and  somewhat  isolated 
radiative  states  with  more  sideways  radiative  transitions. 

In  THz  LO-phonon  designs,  the  approach  is  similar  to  the  first  QCL  developed” 
(which  was  a  mid-IR  laser)  in  that  the  lower  lasing  state  is  depopulated  by  electron- 
phonon  scattering.  Though,  no  terahertz  LO-phonon  QCL  was  fabricated  that  used  all 
diagonal  transitions  like  the  first  mid-IR  QCL.  Also,  the  use  of  digitally  graded  alloy 
injectors  has  been  replaced  with  easier  to  grow  so  called  funnel  or  resonant  tunneling 

'in 

injectors.  ’  THz  designs  differ  in  that  the  radiative  energy  spacing  is  smaller  than  the 
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LO-phonon  energy  levels  spacing.  The  first  THz  QCL  that  used  the  LO-phonon  approach 
for  depopulation,  achieved  the  small  radiative  energy  spacing  by  coupling  the  first  two 
quantum  wells  with  a  thin  barrier.  This  design  approach  still  principally  used  three 
quantum  wells,  but  the  third  well  differed  in  that  it  was  used  to  resonantly  tunnel  the 
lower  lasing  state  (or  equivalently  the  upper  phonon  state)  and  also  arrange  a  lower 
ground  state  spaced  near  the  LO-phonon  energy  (a  fourth  well  was  also  included  at  the 
injector).  Although  some  variations  to  these  structures  have  been  made,  all  of  the  THz 
LO-phonon  based  structures  have  essentially  been  within  the  same  framework  and  have 
included  this  type  of  three  well  arrangement,  differing  only  in  the  design  of  the  injector. 


Figure  I.A.l.  Conduction  band  diagrams  illustrating  the  two  main  broadly  categorized 
QC  structure  architectures,  with  one  section  outlined,  (a)  miniband  structure  (showing  a 
bound  to  continuum  structure),  and  (b)  LO-phonon  structure  (showing  a  three  well 
structure). 


3 


This  three  well  arrangement  (Figure  I. A.  1(b))  is  exemplified  in  the  given  reference. 
Although  at  some  THz  frequencies,  reported  miniband  designs  have  achieved  the  highest 
output  power  levels,  this  is  attributable  likely  only  to  differences  in  the  waveguide 
structure  or  the  number  of  sections,  and  not  the  design  approach  of  the  active  region.  In 
fact,  the  highest  output  power  QCL  used  a  LO-phonon  structure.10 

LO-phonon  based  structures  are  of  interest  because  of  the  large  separation  of  the 
LO-phonon  energy  spacing,  and  also  because  the  highest  output  power  THz  QCL  used  a 
LO-phonon  approach.  This  approach  seemingly  is  more  likely  tolerant  to  thermal 
backfilling,  which  could  be  important  for  increased  operating  temperatures.  The  LO- 
phonon  structures  typically  have  oscillator  strengths  (~  0.5  to  0.96),  while  miniband 
structures  such  as  the  bound  to  continuum  designs  have  higher  listed  oscillator  strengths 
(~  1.9).  However,  this  may  be  somewhat  offset  since  the  LO-phonon  sections  are 
typically  about  half  the  length  of  their  miniband  counterparts,  and  thus  typically  have 
about  twice  the  number  of  (LO-phonon  based)  sections  within  the  same  overall  active 
region  thickness.  Our  multiple  section  analysis  of  some  miniband  structures,  has 
indicated  lower  oscillator  strengths  than  reported.  All  of  these  miniband  and  LO-phonon 
designs  have  used  square  quantum  wells  that  are  symmetric  when  not  under  bias  of  an 
applied  electric  field.  Next,  a  new  theoretical  approach  that  was  investigated  as  part  of 
this  work,  to  the  design  of  the  active  region  for  THz  QCLs,  will  be  discussed. 

B.  STEP  WELL  QC  STRUCTURES 

In  addition  to  the  use  of  miniband  and  square  well  LO-phonon  approaches,  the 
use  of  step  well  (asymmetric  when  unbiased)  QC  structures  was  first  proposed  by 
Freeman  (author)  and  Karunasiri.40  These  types  of  structures  can  allow  for  the  radiative 
and  LO-phonon  transitions  to  be  placed  within  the  same  well.  Since  the  bound  state 
energy  in  high  barrier  square  quantum  wells  increases  ~  n“  (with  n  being  the  quantum 
number  of  the  state),  it  is  not  possible  to  have  an  upper  radiative  THz  energy  spacing  and 
LO-phonon  energy  spacing  (below)  within  a  single  square  quantum  well.  Under  a  step 
well  arrangement,  there  are  three  main  energy  levels,  where  the  transition  from  the  upper 
state  to  the  middle  state  is  at  the  THz  radiative  spacing  and  the  transition  from  the  middle 
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state  to  the  lower  state  is  at  or  near  the  LO-phonon  energy  (~  36  meV  in  GaAs).  The 
middle  state  (upper  phonon  or  lower  lasing  state)  is  a  single  energy  state,  contrasting  to 
previous  LO-phonon  based  QCL  designs  that  have  doublet  states.  Allowing  for  a  single 
middle  energy  state  could  be  important  for  longer  wavelength  (lower  THz  frequency) 
devices  to  reduce  unwanted  absorption.  By  having  vertical  radiative  and  LO-phonon 
transitions  within  the  same  well,  it  is  possible  for  these  types  of  structures  to  yield  high 
oscillator  strengths,  which  can  lead  to  increased  gain  in  the  active  region  provided  the 
upper  state  lifetime  and  injection  efficiency  are  maintained.  The  step  in  the  well  allows 
for  high  injection  efficiency  due  to  the  spatial  separation  of  the  wavefunctions. 

A  step  quantum  well,  in  which  at  least  two  different  conduction  band  heights 
within  a  well  are  used,  allows  for  additional  freedom  and  breaks  the  restriction  single 
square  wells  have.  Thus,  a  THz  and  LO-phonon  energy  spacing  can  be  arranged  within  a 
single  step  quantum  well.  Others  have  analyzed  step  wells  for  proposed  CO2  pumped 
THz  laser  applications.41"42  Those  step  wells  differ  from  the  step  wells  considered  in  this 
research,  where  the  radiative  state  is  positioned  above  the  LO-phonon  transition  states, 
and  are  intended  for  electrically  pumped  QCL  structures.  It  was  these  unique 
characteristics  of  a  step  well,  that  prompted  investigation  of  using  a  step  well  approach 
for  a  QC  structure. 

C.  OVERVIEW 

The  ongoing  imaging  research  using  THz  QCLs  at  NPS  and  by  other  groups,  has 
shown  there  is  a  need  for  high  power  QCL  sources.43"50  The  research  documented  in  this 
dissertation  focused  on  investigating  high  gain  THz  QC  structures  of  a  particular  type, 
LO-phonon  QC  structures.  This  was  because  our  step  well  structures,  which  are  LO- 
phonon  designs,  showed  characteristics  favorable  for  a  high  gain  active  region,  and 
because  the  previously  referenced  highest  power  THz  QCL  was  a  LO-phonon  design  as 
well. 

In  the  sections  that  follow,  techniques  developed  for  analyzing  QC  structures  will 
be  discussed.  Methods  for  modeling  quantum  well  structures,  which  are  an  absolute 
necessity  for  analyzing  or  designing  QC  band  structures,  will  be  covered.  Various 
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scattering  mechanisms  and  the  calculation  of  pertinent  scattering  rates  for  detennining 
state  lifetimes,  will  then  be  discussed.  Fundamental  rate  equations  will  be  also  reviewed 
along  the  way.  As  will  be  shown,  the  scattering  rates  of  the  electronic  states  are  not 
constant,  but  are  density  and  temperature  dependent.  Because  of  this,  rate  equations  are 
not  well  suited  for  analyzing  complete  structures,  as  they  would  be  nonlinear.  One 
approach  that  is  suitable  for  modeling  entire  structures  is  one  based  on  Monte  Carlo 
simulations.  The  Monte  Carlo  method  that  was  developed  for  modeling  the  electron 
transport  in  QC  structures  is  then  discussed.  The  last  sections  cover  the  design  and 
analysis  of  step  well  QC  structures,  and  the  first  sample  that  was  grown  and  processed. 
Comparisons  are  made  to  conventional  square  well  LO-phonon  designs.  The  final  section 
will  show  our  Monte  Carlo  analysis  of  a  square  well  high  power  THz  LO-phonon  QC 
structure,  including  the  effects  of  interface  roughness  scattering,  and  comparisons  are 
made  to  experimental  data. 

The  List  of  Symbols  page  shows  some  of  the  constants  and  notations  used.  Many 
of  the  derivations  work  in  units  with  h  =  c  =  1.  However,  all  sections  dealing  with 
equations  of  final  computable  values  will  have  h  and  c  inserted  back  into  the  equations. 
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II.  QUANTUM  WELLS 


The  energy  states  in  the  active  region  of  a  QCL  are  formed  by  quantum  wells  in 
the  conduction  band.  In  the  sections  that  follow,  methods  for  calculating  the  electronic 
states  in  quantum  well  structures  are  discussed.  The  transfer  matrix  method  is  developed 
for  solving  the  energy  levels  and  wavefimctions  in  such  structures,  both  for  arbitrary 
shaped  potentials  and  for  a  linearly  changing  potential.  The  former  is  useful  for  analyzing 
doped  structures,  where  bending  of  the  band  structure  occurs  due  to  the  separation  of  the 
charges  from  the  host  ions.  This  requires  a  self  consistent  solution  to  Schrodinger’s  and 
Poisson’s  equations.  These  techniques  are  used  to  solve  step  quantum  well  structures,  to 
illustrate  their  potential  use  for  THz  QC  structures. 

A.  TRANSFER  MATRIX  METHOD 

There  are  a  number  of  ways  to  solve  Schrodinger’s  equation  for  bound  states  in 
quantum  well  structures.  One  commonly  used  approach,  the  transfer  matrix  method,  can 
be  used  to  solve  arbitrary  shaped  conduction  band  potentials,  which  is  useful  when  the 
self  consistent  solution  is  needed.  Also  to  be  discussed,  is  the  case  where  the  potential 
can  be  modeled  as  perfectly  linear  (classically  a  charge  moving  under  constant 
acceleration). 

The  most  general  solution  can  be  formulated  by  considering  the  solution  for  the 
case  where  the  potential  is  regionally  constant.  With  this  approach,  any  linear  or 
nonlinear  potential  can  be  modeled  by  a  succession  of  small  steps.  The  solution  to  the 
Schrodinger  equation  I/\i//)  =  E\  y/)  for  a  region  with  constant  potential  can  be  written  as 
a  sum  of  forward  and  reverse  waves  (Figure  II.A.l). 

y/{x)  =  Aeikx  +  BAQikx  (II.A.l) 

Since  the  wavefuncion  must  be  continuous,  at  an  interface  (x„)  the  first  boundary 
condition  can  be  written  as 

¥nixn)  =  Vn+ i(*J  (II.A.2) 
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Figure  II.  A.  1 .  Illustration  representing  the  conduction  band  as  a  number  of  constant  steps, 
showing  the  coefficients  for  use  in  the  transfer  matrix  method. 


A  second  boundary  condition,  is  that  the  first  derivative  must  be  continuous.  Choosing  to 
take  into  account  the  change  in  effective  masses  across  the  boundaries,  it  is  noted  that  the 
conserved  probability  density  flux  is  j  =  Mm  Im{  i//V  i//} .  Setting  /,,  =  /'„+/  the  following  is 
found 


1  dyfn{xn)  _  1  ch//n+l(xn) 

mn  dx  m„+1  dx 


(II.A.3) 


The  transfer  matrix  can  be  found  from  these  boundary  conditions  by  seeking  relations  for 
the  coefficients  A„  and  Bn  in  terms  of  An+ 1  and  Bn+\.  Substituting  equation  (II. A.  1)  into 
equations  (II.A.2)  and  (II.A.3),  the  relationships  for  An  and  Bn  can  be  found.  The  resulting 
equations  can  be  cast  into  matrix  fonn 
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where  inserting  fr  back  into  the  equations  k„  =  [2mn(E-Vn)/h  ]  .  Because  of  the  way  the 
transfer  matrix  is  developed,  correct  boundary  conditions  at  the  interfaces  are  always 
insured.  In  a  structure  with  N  regions,  the  matrix  equation  in  terms  of  the  first  and  last 
regions  may  be  found  as 


1 

o 

_ 1 

ii 

s 

o 

s 

1 

to 

~An- ll 

M 

Vil 

1 

3 

o 

o 

moi" 

- 1 

_ 1 

LV 

A-lJ 

_mio 

mll  _ 

(II.A.5) 


For  bound  state  solutions,  the  wavefunctions  must  be  square  integrable  and  therefore 
vanish  as  x  — >  ±oo.  From  equation  (II.A.5)  we  see  the  bound  eigenstates  can  be  found  by 
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setting  mooCS)  =  0.  The  transmission  coefficient  can  also  be  found  for  the  case  where  an 
electron  is  tunneling  through  one  or  more  barriers.  Assuming  an  electron  incident  from 
the  left  means  that  Bn- i  =  0.  The  transmission  coefficient  is  then  found  from  T(E)  =  y ,v- 1  //o 
=  kN-imo/(komN-i)\l/moo\  . 

If  the  structure  can  be  modeled  as  having  a  perfectly  linear  potential  (as  is  the 
case  for  some  QC  structures  where  minimal  bending  of  the  band  structure  occurs  from 
the  doping),  the  transfer  matrix  can  be  formulated  for  the  case  of  a  linear  potential. 

Writing  the  potential  as  V(x)  =  Vn-e\E\x,  Schrodinger’s  equation  becomes  Airy’s  equation 

2 

2  =pn(x)vn(x)  (II.A.6) 

dpn  (x) 

2  2  1/3 

where  pn(x)  =  (— |E| x-r/n)/j3,  ijH  =  (E-Vn)/e,  ft=  [h~\E\  /(2 m,,e)\  ,  and  the  solutions  are  a 
linear  combination  of  Airy  functions. 


Vn  0)  =  AnAi(Pn  O))  +  BnBi(Pn  O))  (II.A.7) 

The  transfer  matrix  can  be  obtained  in  a  similar  manner  as  with  the  constant  potential 
case,  and  can  be  simplified  to51'52 


A 

m2n 

Ai(an+l)Bi  '(a,,)  "2/3  Ai  '(an+1)Bi(an) 

m213 

Bi(an+1)Bi  '(«„)  ”/3  Bi(a„)Bi  '(«„+1) 

An+ 1 

™n+X 

>K+i 

=  n 

2/3 

2/3 

”2/3  '(an+ 1)  Ai  '(an)Ai(an+l) 

l>K+i 

i'll 

"2/3  Ai(an)Bi  '(a„+1)  Ai  '(a„)Bi(an+l ) 

_Bn+ 1. 

(II.A.8) 

where  a„  =  (2m„e\E\/h  )  (—x„— )]rJ\E\)  and  the  prime  denotes  a  derivative.  In  a  structure 
with  N  regions,  we  can  again  write  the  matrix  equation  in  terms  of  the  first  and  last  region 
the  same  as  in  equation  (II. A. 5).  For  bound  state  solutions,  the  wavefunctions  must  again 
vanish  as  x  — >  ±oo,  or  strictly  speaking  prior  to  the  classical  turning  point  as  x  tends 
toward  oo.  Since  Bi(an )  — >  oo  as  x  — >  -oo,  we  set  B0  =  0,  noting  that  mi0  ~  0  and  mi0  « 
mu.  The  bound  eigenstates  are  found  by  setting  m n(E)  =  0. 

1.  Self  Consistent  Solution  of  Schrodinger’s  and  Poisson’s  Equations 
The  self  consistent  solution  of  Schrodinger’s  and  Poisson’s  equations  is  necessary 
to  take  into  account  doping  in  quantum  well  nanostructures.  Formulating  the  transfer 
matrix  method  for  a  linear  bias  is  computationally  efficient,  but  is  technically  only  valid 


9 


for  undoped  structures.  In  order  to  take  into  account  bending  of  the  band  structure  that 
occurs  due  to  the  separation  of  the  dopant  charges  from  the  host  ions,  Schrodinger’s 
equation  is  solved  by  approximating  the  conduction  band  as  a  number  of  small  steps 
using  the  transfer  matrix  method.  This  approach  allows  for  solution  of  any  arbitrarily 
shaped  potential.  Another  popular  alternative  technique  used  is  the  shooting  method, 
which  deals  with  the  effective  mass  at  intermediate  points  by  taking  the  mean  at  the 

53 

neighboring  points. 

The  solution  to  Poisson’s  equation  can  be  fonnulated  for  numerical  computation 
as  follows.  Starting  from  Gauss’s  law  V-E  =  pi  8  and  the  electrostatic  form  of  Faraday’s 
law  VxE  =  0,  noting  that  the  curl  of  the  gradient  always  vanishes,  the  electric  field  can  be 
written  in  the  usual  way  as  E  =  -V Va.  This  leads  to  Poisson’s  equation  V"Va=  -pis.  Now 
consider  in  general  the  following  infinite  sheet  shown  in  Figure  II.A.2,  which  is  a  plane 
surface  of  thickness  Ax3  with  charge  density  mo(x),  illustrating  the  electric  field. 


side  2 


side  1 


Figure  II.A.2.  Infinite  plane  surface  of  thickness  Ax3  and  charge  density  «3d(x),  showing 
the  electric  field. 


Applying  Gauss’s  law  we  find 
f d2xE  ■  n  =  —  [ d3xp(x)  =—  f  d2xa(x) 

J  r»  J  c  J 


(II.A.9) 


where  cr  is  the  2D  charge  density.  For  this  case,  the  electric  field  is  then  found  from 
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(II.A.10) 


a 
s 

The  charge  density  as  a  function  of  the  growth  direction  may  be  written  as 
cr(x3)  =  -e[n2D^'*(x3)^r(x3)  -  «3d(x3)]Ax3  (II.A.11) 

Using  equations  (II.A.10)  and  (II.A.11),  the  electric  field  is  found,  and  the  potential  can 
then  be  determined  by  integration 

Vr7=-jdx  En  (II.A.12) 

If  the  dopants  are  concentrated  in  other  states  than  just  the  ground  state,  the  following 
modification  is  needed 

a(x3)  =  -e£  [nmiy/*  (x3)<//.(x3)  -  «3d,/(u)Ku  (II.A.13) 

i=0 

where  the  sum  of  all  the  sheet  densities  for  all  of  the  states  equals  the  total  sheet  density. 
The  self  consistent  solution  iteratively  computes  the  perturbing  potential  and  adds  it  to 
the  conduction  band  profile  until  after  a  given  number  of  iterations,  a  self  consistent 
solution  is  reached.  Equation  (II.A.13)  assumes  the  charge  distributions  in  the  states  are 
known.  In  practice,  this  is  not  known  unless  the  electron  transport  of  the  full  QC  structure 
is  analyzed  and  solved.  One  method  for  solving  the  state  populations  can  be 
accomplished  using  Monte  Carlo  simulations,  and  this  will  be  discussed  later. 

For  QC  structures,  doping  is  needed  to  introduce  charge  carriers  within  the  active 
region,  and  also  to  help  reduce  the  possibility  of  the  formation  of  high  field  domains.54  A 
typical  LO-phonon  QC  structure  may  have  on  the  order  of  ~  200  cascaded  sections,  and 
an  undoped  active  region  would  essentially  represent  a  large  intrinsic  region.  Minimal 
doping  levels  can  result  in  lower  threshold  currents,  and  moderate  doping  levels  have 
been  shown  to  result  in  the  highest  operating  temperatures.55,56  In  practice,  only  one  well 
or  a  portion  of  a  well  or  barrier  is  doped,  in  each  section  of  the  cascade. 

B.  STEP  QUANTUM  WELLS 

As  previously  discussed,  a  step  quantum  well  is  not  limited  to  have  the  bound 
state  energy  increase  ~  n“  as  is  the  case  for  a  conventional  square  quantum  well.  As  such, 
the  THz  and  LO-phonon  energy  spacing  can  be  arranged  within  the  same  step  well.  This 
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approach  is  illustrated  using  the  step  quantum  well  structures  shown  in  Figure  II.B.l. 
These  wells  are  comprised  of  AlxGai_xAs  layers  with  compositions  of 
0.143/0.035/0/0.143.  The  step  well  thicknesses  in  mn  are  29.9/12.3  for  Figure  II. B.  1(a) 
and  20.9/13.5  for  Figure  II.B.l(b).  The  energy  difference  between  states  2  and  1  for 
Figures  II. B. 1(a)  and  (b)  are  £21  =  13.7  meV  (3.3  THz)  and  £21  =  17.9  meV  (4.3  THz) 
respectively.  Both  step  wells  have  El0  near  the  LO-phonon  energy.  The  energy  levels  and 
wavefunctions  were  calculated  assuming  an  applied  field  of  10.1  kV/cm  (typical  for  LO- 
phonon  QCLs)  using  the  material  parameters  given  in  the  references.57’65  This  illustrates 
that  it  is  possible  to  arrange  different  radiative  transition  energies  and  keep  near  resonant 
phonon  transitions,  within  a  single  step  well.  This  approach  can  be  used  to  design 
radiative  frequencies  lower  and  higher  than  these  as  well. 


Figure  II.B.l.  Conduction  band  profiles  of  step  quantum  wells,  (a)  Conduction  band 
profile  of  a  step  quantum  well  comprised  of  AlxGai  xAs  layers  with  compositions  of 
0.143/0.035/0/0.143  and  well  thicknesses  in  nm  of  29.9/12.3,  with  £21  =  13.7  meV  (3.3 
THz)  and  £10  =  37.9  meV.  (b)  Conduction  band  profile  of  a  step  quantum  well  comprised 
of  AlxGai_xAs  layers  with  compositions  of  0.143/0.035/0/0.143  and  well  thicknesses  in 
nm  of  20.9/13.5,  with  £21  =  17.9  meV  (4.3  THz)  and  £10  =  36.5  meV.  For  both  (a)  and 
(b),  the  1  to  0  transitions  are  near  the  LO-phonon  energy  and  the  applied  bias  is  10.1 
kV/cm. 


This  approach  has  a  number  of  attractive  features.  Since  the  radiative  and  LO- 
phonon  transitions  are  both  vertical,  this  approach  can  yield  large  oscillator  strengths  and 
fast  LO-phonon  scattering  rates  for  depopulation.  Oscillator  strengths  for  step  well  QC 
structures  are  typically  around  unity  or  greater.  Conventional  square  well  LO-phonon 
structures  typically  have  oscillator  strengths  only  about  half  of  this  and  up  to  0.96.  It  is 
also  noted  that  the  step  represents  an  additional  barrier  that  can  reduce  parasitic  injection 


12 


into  the  middle  state.  This  could  help  the  injection  efficiency  in  QC  structures,  especially 
for  longer  wavelength  devices  where  the  radiative  energy  transition  spacing  is  small.  It  is 
also  seen  that  the  middle  state  is  a  single  state,  contrasting  to  other  conventional  square 
well  THz  LO-phonon  designs  that  always  have  doublet  middle  states.  Having  a  single 
state  could  help  reduce  unwanted  THz  absorption,  which  could  be  important  for  longer 
wavelength  devices.  In  the  following  sections,  the  calculation  of  various  scattering  rates 
of  the  electronic  states  in  quantum  well  structures  will  be  discussed. 
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III.  SCATTERING  AND  OPTICAL  TRANSITION 


In  order  to  determine  the  state  lifetimes  within  QC  structures,  various  scattering 
rates  of  those  states  must  be  computed.  In  the  sections  that  follow,  the  calculation  of 
scattering  rates  for  the  pertinent  scattering  mechanisms  will  be  covered.  As  will  be 
shown,  the  scattering  rates  and  hence  the  corresponding  lifetimes  of  the  electronic  states 
are  not  constant,  but  rather  are  density  and  temperature  dependent.  We  will  also  discuss 
the  optical  transition,  and  the  important  calculation  of  the  gain.  Some  approximate  rate 
equations  will  also  be  used  to  illustrate  fundamental  device  parameters  of  interest. 


A.  ELECTRON-PHONON  SCATTERING 


Phonons  are  quantized  vibrations  within  a  lattice.  Many  electrical  and  thermal 
properties  of  materials  are  governed  by  phonon  interactions.  In  semiconductor 
heterostructures,  electron-phonon  scattering  can  often  be  the  dominant  scattering 
mechanism.  If  the  spacing  between  two  states  is  at  or  near  the  resonant  LO-phonon 
energy,  the  scattering  rate  to  the  lower  state  can  be  very  fast.  This  is  why  resonant 
phonon  scattering  has  been  used  in  LO-phonon  QCL  devices  as  the  primary  mechanism 
for  keeping  a  population  inversion.  Intraband  LO-phonon  scattering  along  with  electron- 
electron  scattering,  are  the  mechanisms  that  cool  the  electron  gas  within  subbands  and 
thermally  distribute  the  electrons  into  Fermi-Dirac  distributions.  Since  our  focus  is  on 
LO-phonon  QC  structures,  this  is  the  most  important  scattering  mechanism,  and  it  will  be 
discussed  first. 

Beginning  by  considering  a  simple  classical  diatomic  system,  which  is  a  good 
model  for  GaAs,  the  model  consists  of  masses  attached  via  springs  in  a  chain  (Figure 
III.A.  I).66  The  Hamiltonian  for  the  system  is 

h  =  Yjke  +  v 


I 


P 2  P: 

X  :  X  : 

+  ■ 


2ml  2m2 


(III.A.  1) 
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Figure  III.A.l.  Linear  chain  of  diatomic  atoms  of  masses  m\  and  m2. 


From  Hamilton’s  equations  we  can  write 
dH 


P  =  — 


8x 


mxXj  =  - 

m2x'j  =  - 


k{x 


—  K\Xj  —X. 


k{. 


M 

)-*(• 


Xj-Xj 


X:-X1)-K\Xj^-X; 


;)] 


(III.A.2) 


or 

mxXj 

=  k 

+x'j_  1  - 

2xj) 

m2x'j 

=  } 

'c(xj+l+Xj- 

-2  x'j) 

We  seek 

:  solutions  of  the  form 

xj  = 

Re 

OJ 

11 

X 

Re 

\X'  Q^,)a~<a^ 

} 

(III.A.3) 


(III.A.4) 


where  a  is  the  lattice  spacing  between  like  atoms.  Using  equations  (III.A.4),  equations 
(III.A.3)  can  be  written  in  matrix  notation  after  some  manipulation  as  follows  (Re{}  is 
assumed) 


co2 

(  \ 

~X~ 

(nr  —  co~) 

m\ 

-(l  +  e  ,qa ) 

\  m  1  / 

2 

\  / 

2 

-  (°m2  (l  +  e'qa ) 

2  v  ’ 

[py2 

\  ml 

-co2) 

X' 

where 


(III. A. 5) 


2  2  k 

1  =  - 

mi  mx 

2  2k 


(III. A. 6) 
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and  q  =  qrqr\.  The  solution  is  found  by  setting  the  determinant  to  zero.  Solving  the 

2 

quadratic  in  terms  of  co  leads  to  the  following. 


co1  = 


(co2  +  co1 

\  mx  m 


-  + 
1 


(  2  ,  2  \2 

0)m  +  Q) 


2  2-2 

-  (0m  0)m  sin 


^  qa^ 


V  z  J 


(III. A. 7) 


Equation  (III. A. 7)  can  be  used  to  graph  the  dispersion  curve  of  co  versus  q.  Shown  in 
Figure  III.A.2  is  the  dispersion  curve  with  m2  =  0.93  in  i  (mGa  =  69.72  g/mol,  /«as  =  74.92 
g/mol).  As  can  be  seen,  there  are  two  main  branches  of  phonons,  known  as  acoustic  and 
optical  phonons  (the  notation  co  — >  coAc  for  the  acoustic  branch  and  co  — *  coop  for  the 
optical  branch,  will  be  used). 


Figure  III.A.2.  Phonon  dispersion  graph,  showing  the  acoustic  and  optical  branches. 


For  the  acoustic  branch,  which  corresponds  to  the  “negative  sign”  root  in  equation 

1/2 

(III. A. 7),  as  q  — >  0,  co  — >■  0  and  we  see  using  sin(go/2)  ~  qa/2  and  (1  +  stuff)  ~  1/2 (stuff) 
qa 


co  ac 


co  co 

ml  m2 


(III. A. 8) 


Substituting  equation  (III.A. 8)  into  equation  (III.A. 5)  gives  the  relative  displacement 
between  the  two  atoms. 


Re 


X\  =  X' 


(III. A. 9) 
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For  the  optical  branch,  which  corresponds  to  the  “positive  sign”  root  in  equation 
(III. A. 7),  as  q  — >  0,  co  — >  nw,  and  when  q  =  0  the  following  is  found 


(III.  A.  10) 


where  l/mredUcded  =  l/mi+l/»i2.  Substituting  equation  (III.A.10)  into  equation  (III.A.5) 
gives  the  relative  displacement  between  the  two  atoms. 


(III.  A.  11) 


As  expected  because  of  the  high  frequency  nature  at  q  =  0,  the  optical  branch  is 
characterized  by  opposing  motion  of  the  opposite  atoms. 

From  equation  (III.  A. 7)  it  is  seen  that  there  are  a  range  of  solutions  specified  by  q 
=  ±n/a.  For  periodic  boundary  conditions  (x,+n  =  Xj  and  x]+n  =  xj)  of  a  finite  chain  of 
length  L  =  Na,  i.e.,  the  ends  of  the  chain  are  tied  together  in  a  ring,  qL  =  2 Tin,  or  q„  = 
27vaJL.  There  are  N  number  of  ni\  atoms  =  N  number  of  m2  atoms  =  N  number  of  two- 
body  pairs  or  2 N  number  of  atoms  total  (which  is  the  number  of  degrees  of  freedom  for 
the  system).  For  the  3D  case,  in  addition  to  longitudinal  phonons,  there  can  also  be 
transverse  phonons.  In  mixed  composition  crystals,  the  main  branches  can  split  into  sub¬ 
branches.  This  has  been  experimentally  measured  and  theoretically  derived  using  the 
Lyddane-Sachs-Teller  (LST)  splitting.67'69 

In  polar  semiconductors,  the  phonons  create  polarizations  and  thus  electrostatic 
fields  that  follow  the  phonon.  The  interaction  that  describes  this  is  referred  to  as  the 
Frohlich  interaction.  An  expression  for  this  interaction  can  be  found  by  first  considering 
in  general  the  dipole  moment  pi  =  qc\,  where  qc  is  the  charge  and  |x|  is  the  displacement. 
The  dipole  moment  per  unit  volume  V,  referred  to  as  the  electric  polarization  P,  is  then 
written  as 


P  =  Iy(Pi)  (III.  A.  12) 

The  induced  polarization  due  to  the  deformation  by  a  phonon  can  be  written  as  (noting  e 
is  the  deformation  charge,  i.e.,  the  charge  that  is  transferred,  not  the  electron  charge) 
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(III.  A.  13) 


P  =e*(X-X')— - 
V  V 

The  continuity  equation  djv=  0,  where  J' '  =  (p,  J),  can  be  integrated  with  respect  to  time 
to  give 


J  dtV  J  =  -^dt- 


(III.  A.  14) 


Since  P  =  \dtJ,  the  charge  density  (charge  per  unit  volume)  is  p  =  -VP.  For  an  optical 
phonon  with  wavevector  q,  this  becomes 

p(\,  t)  =  Rej-  eV  ■  (X  -  X')—  =  Re{-  f eiq  ■  (X  -  X')— (III.A.  15) 


The  dot  product  q-(X-X')  projects  only  in  the  longitudinal  direction,  i.e.,  the  q  direction, 
thus  only  longitudinal  optical  (LO)-phonons  can  induce  a  charge.  The  induced 
electrostatic  potential  (O^)  can  now  be  found  from  Poisson’s  equation 
V“® q  =  -pis,  where  s=  and  noting  from  the  divergence  theorem 


J</2xV2®?  =-|j 


=  -[d2x-?~ 


(III.A.16) 


O,  =  ReJ  f  (Z  -  X%  eud2u  \  -  Re]  f ‘Px  -  X'Pj 


V  (iq) 


_^7V( 
£\  V 


(III.A.  17) 


where  u  =  /qx  -ia>qt,  du  =  iq-dx,  and  d2u  =  /q  d2x.  Perfonning  the  integral  and  taking 
the  real  part,71  the  following  is  found 


®  — -(A-X')(e'qx  -e~'qx)  =  ————{X -X')(-e,qx  +  e~'qx)  (III.A.  18) 

2  en  V  q  2  en  V  q 

and  the  perturbing  potential  is  now  found  from 


e  -phonon 


q  2  v  q 


(III.A.  19) 


-  e'q  x  +  e~'q : 


An  expression  for  the  deformation  charge  e  due  to  a  LO-phonon  can  be  determined  from 


classical  arguments,  and  is  found  to  be6 
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e  = 


I  mV  0JW 


N 


1  1 


V^O 


(III.A.20) 


'st  J 


To  solve  for  the  amplitude  X-X\  the  system  must  be  quantized.  Second 
quantization  could  be  used  to  promote  everything  directly  in  our  Hamiltonian  to  operators 
using  the  canonical  commutation  relations,  but  there  is  an  easier  way.  This  standard 
procedure  is  to  write  the  Hamiltonian  in  terms  of  generalized  coordinates.66  Using 
canonical  coordinates  and  momenta  (Q/  and  P / ),  a  set  of  normal  coordinates  (qq)  is  found 
such  that  the  potential  is  in  diagonal  fonn.  The  canonical  commutation  relations  can  then 
be  used  to  quantize  the  system  in  a  straight  forward  fashion.  The  Hamiltonian  expressed 
in  canonical  coordinates  is 


"*»=Zt-+ZtVW; 


'2  rn, 


where  for  the  restricted  harmonic  oscillator  case 

82V 


A°  SQ'SQ' 


rQ'Qj 


(III.A.21) 


(III.A.22) 


noting  that  Ay  is  real  and  symmetric  (Ay  =  A /;)  and  [Q/,P/]  =  i5y,  [Q/,Q/]  =  [P/,  P/]  =  0. 
Following  Feynman72  by  scaling  our  canonical  coordinates  and  momenta  (to  make  the 
notation  simpler),  we  look  to  find  the  normal  coordinates.  The  scaling  chosen  is  Qj  = 
(m,)12 Q/,  Pj=  P *j/(mi)m  and  Uy  =  2A y  /(m,mj)12.  The  Hamiltonian  then  becomes 

Hm  =  \'LPA\YU„Q,Q,  (III.A.23) 

^  i  ^  ij 

noting  that  Uy  is  real  and  symmetric  (Uy  =  Uji).  A  transform  between  the  nonnal 
coordinates  and  the  canonical  coordinates  can  be  written  as 

%  =  Zc,/0  (III.A.24) 

i 

where  the  transformation  matrix  Cq(-  is  orthogonal,  ^CaiCfJi  =  Sq/3,  ^Cqt.Ca/  =  8ih  and  the 

i  q 

inverse  transform  is 

0=  Z0A  (III.A.25) 

q 

The  canonical  momenta  is  likewise  written  as  (where  \qq,pp[  =  iSqp) 
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(III.A.26) 


p,=Zc^ 


The  eigenvalues  o\  can  now  be  found  from 


and  thus 


nuijQiQj=i:« 


(III.A.27) 


(III.A.28) 


Using  equations  (III.A.26)  and  (III.A.28),  the  Hamiltonian  can  now  be  written  as  a  sum 
of  independent  harmonic  oscillators. 


=  +  %% 


Z  q  Z 

(  2  .  22^ 


(III.A.29) 


These  decoupled  oscillators  represent  noninteracting  phonons.  The  solution  is  of  course 
well  know  for  the  harmonic  oscillator,  and  it  is  now  straight  forward  to  quantize.  Treating 
each  classical  hannonic  oscillator  as  a  quantum  harmonic  oscillator,  pq  and  qq  are 
promoted  to  operators,  and  the  canonical  commutation  relations  are  imposed. 
q  ,  n  =  id 

1  r  "  1  (III. A. 30) 

<ii,qj\=[Pi,Pj\=0 


The  operators  p  and  q  can  be  written  in  terms  of  ladder  operators,  noting  that  the 
canonical  commutation  relation  [q,,p,]  =  iSy  is  the  same  as  [aq,aq]  =  1. 


(III.A.31) 


The  Hamiltonian  for  the  sum  of  the  independent  oscillators  can  now  be  written  as 


HHO  =Z^q[aiaq+|]  =  Z^[aqa!“|j  (III. A. 32) 

and  the  Hamiltonian  can  be  explicitly  written  for  the  3D  case  as  follows,  keeping  in  mind 
that  for  the  3D  case  both  longitudinal  and  transverse  modes  can  exist 
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3  f  l 

^0=ZZ®J  aLaq,.  +  o 


q  5=1 


(III. A. 33) 


where  for  each  independent  oscillator  |i  =  0)  is  an  eigenstate  of  Hho  with  a  zero-point 
eigenvalue  of  Am,,.  When  the  zero-point  energy  is  set  to  zero,  aqjS|i  =  0)  =  0  and  the 
vacuum  state  |0)  has  an  energy  eigenvalue  of  £  =  0.  Using  the  commutators  below 


HHO,  aq,.s] 


^qaq,.s 


[Hm,  aq,.s] 


^qaq  ,.s 


(III. A. 34) 


it  is  seen  |n)  =  (aqjS)n|0)  and  the  Hamiltonian  for  each  oscillator  is  diagonalized,  which  has 
a  complete  set  of  eigenvalues  of  (n+l/2  )a>=  (n +V2)ELo  (for  n  =  0,  1,  ...). 

The  background  energy  of  the  phonon  modes  is  224  <uq.  For  systems  with  many 
degrees  of  freedom,  this  c-number  becomes  large  and  infinite  in  extent  for  fields.  Since 
only  energy  differences  can  physically  be  measured,  this  zero-point  value  is  not 
important  for  our  purposes  here.  However,  it  is  noted  though  for  fields  that  have  infinite 

n'i 

degrees  of  freedom,  this  shift  in  the  zero-point  could  potentially  be  a  problem. 

The  statistics  of  phonon  particles  can  be  determined  by  considering  the  following 
two  particle  state  apaq|0).  Noting  that  ap  and  aq  commute,  the  state  is  identical  if  the  two 
particles  are  interchanged.  Analogous  to  a  classical  oscillator  that  can  be  excited  to  an 
arbitrary  number  of  high  levels,  mode  q  can  have  an  arbitrary  number  of  particles.  Thus, 
phonons  are  Bosons  and  the  average  number  of  phonons  in  a  mode  can  be  found  from 
Bose-Einstein  statistics  ( N)  =  1/Z  Tr  Ac  /,7/. 

The  mode  displacement  is  now  found  by  using  the  correspondence  principle 
where  the  classical  energy  in  the  mode  is  the  same  as  the  quantum  energy  in  the  mode, 
which  is  NqEw  (where  Nq  is  the  average  number  of  phonons  in  the  mode  and  ELo  is  the 
quantized  eigenstate  energy). 


\nuow2 \X  - Xf  N  =  NqEw  (III.A.35) 

\X  -  X'\  =  I2N^Elo  (III. A. 36) 

1  1  y  Nmcow2 

Substituting  the  expressions  for  e*  from  equation  (III.A.20)  and  the  modal 
amplitude  X-X'  from  equation  (III.A.36)  into  equation  (III.A.19),  the  perturbing  potential 


for  the  LO-phonon  interaction  can  be  written  as 
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where  q  =  |q|,  which  is  the  magnitude  of  the  3D  momentum  vector  (includes  the  parallel 
in-plane  and  X3  components).  Separating  this  into  emission  and  absorption  terms,  with  Nq 
—>  Nq+\  for  the  emission  case,  to  include  spontaneous  and  stimulated  phonon  emission 


j.  1  /o 

This  can  also  be  written  in  terms  of  ladder  operators  (a/V/  =  (Nq+ 1 )  |A^+1),  &q\Nq)  = 
{Nq)m\Nq-\)) 


Sometimes  this  is  expressed  in  terms  of  the  nondimensional  constant  otp  = 
e2/(8  nh  co)(2m  colh)m(\l  1  /  sst ) . 

The  following  diagram  (Figure  III.A.3)  shows  the  indexing  and  notation  used  for 
electron-phonon  scattering,  where  i  is  the  initial  electron  state, /is  the  final  electron  state, 
and  n  is  the  subband  index. 
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Figure  III. A. 3.  Feynman  diagram  for  electron-phonon  scattering. 


In  this  section  some  vectors  have  three  components  (q)  whereas  others  have  two  parallel 
in-plane  components.  Parallel  subscripts  will  be  used  for  parallel  in-plane  vectors  (k[|  and 
qn)  and  the  prime  notation  will  indicate  the  final  state.  Now  with  the  Frohlich  interaction 
potential,  the  scattering  rate  can  be  computed.  The  2D  wavefunctions  are  Bloch 
functions. 

|/)  =  \uiai'j  (III.A.40) 


Since  the  interest  is  in  the  matrix  element  for  intersubband  transitions,  the  following 
approximation  can  be  used74 


(vf  \v\  Vi  )  =  {uf  af  \v\ 11  >  ai  )~(uf  \v\ 11  i ) 

f\V  \ai) 


cell  X*'  f 


a  f  a,- )  +  ( u  f  \u 


■f  r< 


cell 


a  f  V  a, 


(III.A.41) 


where  for  intersubband  transitions  (a) a)  =  0  and  (uJ[U,)ceu  =  1.  Working  in  the  (x| 
representation,  where  {x\,  a'2.  x_i)  =  (x,  y,  z)  with  x\,  xi  parallel  to  the  growth  direction  X3, 
the  wavefunctions  can  be  written  as 


x  7  = 


u 


eik|rV„(^): 


VI 


e'k  V„(^) 


(III.A.42) 


This  will  be  used  in  subsequent  sections  and  the  periodic  n,  part  of  the  Bloch  function 
will  be  dropped. 

1.  LO-Phonon  Emission 

Considering  the  emission  case  first,  where  V  =  V  ,  .  The  matrix  element  can 

^  e  -phonon 

be  written  as 
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(III.A.43) 

For  mode  q,  there  is  only  one  value  of  qn  that  makes  the  d  x\  integral  nonvanishing,  and  is 
such  when  the  exponent  is  zero  (otherwise  the  other  values  lead  to  the  addition  of  random 
phase  terms  which  add  to  zero).  This  gives  a  Kronecker  delta  function,  which  is  the 
conservation  of  momentum  as  expected  for  the  emission  case. 
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(III.A.44) 

To  find  the  electron-phonon  scattering  rate,  the  matrix  element  is  substituted  into  Fermi’s 
golden  rule,  summing  over  the  final  states. 7:1  Inserting  h  and  c  back  into  the  equations 
(noting  El0  =  hcoLo ) 


2  n 


w_  ,  =— y  (/|Fii) \s\ 

e  -phonon  ^  I  I  / 1  | 
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(III.A.45) 

Converting  the  final  states  summation  to  an  integral  in  k-space75  and  dropping  the 
momentum  delta  function 
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•  iq  ii 


Figure  III.A.4.  Contour  illustrated  for  closing  around  the  lower  pole,  for  evaluating  the 
contour  integral  in  equation  (III.A.49). 


26 


J  dq 
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- - =  —  e for  (x3  -x^)  >  0  (III. A. 50) 


-iq^+iq^  ql{ 

and  similarly  for  negative  X3-X3,  the  contour  is  closed  in  the  upper  half  plane 


j  dy 


C  TC  -q..  |x3  -x3 1  r  /  t  \  r\ 
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(III. A. 51) 


J  dq 


=  — e 
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(III. A. 52) 


The  /2D(qn)  integral  can  now  be  written  as 

72d  fay )  =  —  f  dx3dx'3i//*n,  (x3  )(//„  (x3  )y/n,  (x'3  )y/n  {x[  )e 
?ii 
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(III. A. 53) 


If  we  wish  to  include  a  screening  term,  equation  (III.A.49)  with  screening  becomes 


76 
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and  the  contour  integral  becomes 
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(III. A. 55) 


-3/9II  +9»c  1*3  “*3  | 


Finally,  upon  substitution  the  /2D(q||,fe)  integral  becomes 
^2D  fa||,^c)  =  dx,dx,\j/*n,  ( x3  )y/n  (x,)i//*t(  x'3  )y/n  ( x( ) 
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(III. A. 56) 


Now  consider  the  energy  delta  function  in  equation  (III.A.47),  which  conserves 
energy.  The  inner  function  of  the  delta  function  can  be  simplified  as  follows 
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(III. A. 57) 


En,+Ek,+h(DLO-En-Ek=  0 


nVw-{E„-En')+Ek’-Ek  =0 

n°>'e  +ek-%  ~Ek  =  0  where  ha>'e  =ficoLO  ~(En  -En) 


(k\\  ~  2k\H\\  C0S  9  +  ?||  )-  Ek  =  0 


where  a  similar  notation  to  the  given  reference  has  been  used.77  Note  that  $  is  the  angle 

between  k||  and  q||.  The  expression  may  be  rearranged  in  tenns  of  the  quadratic. 

ti~  fi~ 
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The  roots  to  the  quadratic  are 


(III. A. 58) 
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2  m  ol. 


(III. A. 59) 


and  the  following  is  noted 
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Using  S{f(x))  =  YJ5{x-xi)/\  f  (xt )  | ,  where  Xi  are  the  roots  of  fix),  the  delta  function 

i 

becomes 
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The  final  expression  for  the  emission  case  is 


(III.A.61) 
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(III.A.62) 


Now  consider  the  different  emission  cases  that  are  possible.  There  are  three 
different  emission  cases  that  can  arise,  and  are  referred  to  as  cases  la,  lb,  and  2. 

LO-phonon  emission  case  la:  Ek  <  hco'e  =  hcow-{En-  En)  and  hco'e  >  0,  i.e.,  E„- 
En'  <  hcow  and  the  electron  does  not  have  sufficient  in-plane  kinetic  energy,  the  transition 
is  forbidden.  The  subbands  are  spaced  such  that  the  energy  difference  is  smaller  than  the 
LO-phonon  energy  and  the  electron  does  not  have  enough  in-plane  kinetic  energy,  so  the 
transition  is  forbidden  (Figure  III. A. 5). 


Picolo 


Figure  III.A. 5.  LO-phonon  emission  case  la,  En-E„’  <  hcoLO  and  the  electron  does  not 
have  sufficient  in-plane  kinetic  energy,  forbidden  transition. 
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LO-phonon  emission  case  lb:  Ek  >  hco'e  =  ha>Lo-(En-En)  and  co'e  >  0,  i.e.,  En-En' 
<  ha>w  and  the  electron  has  sufficient  in-plane  kinetic  energy,  the  transition  is  allowed. 
The  subbands  are  spaced  such  that  the  energy  difference  is  smaller  than  the  LO-phonon 
energy  and  the  electron  has  enough  in-plane  kinetic  energy,  the  transition  is  allowed 
assuming  no  state  blocking  (Figure  III.A. 6).  The  red  and  blue  circles  on  the  E-k  diagram 
are  equal  energy  planes  and  show  the  possible  range  for  the  in-plane  momentum  vectors. 
The  final  state  represented  by  the  blue  circle  must  be  vacant,  otherwise  state  blocking 
will  prevent  such  a  transition.  There  are  two  possible  solutions  for  this  case,  q\\±e.  Limits 
on  0  are  necessary  to  keep  the  momentum  vectors  real  valued  (to  keep  from  getting  a 
negative  sign  under  the  square  root).  Also  note  that  at  resonance  the  transition  is  allowed. 


Figure  III.A. 6.  LO-phonon  emission  case  lb,  En-En'  <  ficow  and  the  electron  has 
sufficient  in-plane  kinetic  energy,  allowed  transition. 

Equation  (III. A. 62)  becomes 
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(III. A. 64) 


LO-phonon  emission  case  2:  Any  Ek  is  valid  and  co'e  <  0,  i.e.,  En-En'  >  Hcolo,  the 
transition  is  allowed.  The  subbands  are  spaced  such  that  the  energy  difference  is  larger 
than  the  LO-phonon  energy,  and  the  transition  is  always  allowed  assuming  no  state 
blocking  (Figure  III.A. 7).  This  “sideways”  transition  has  one  solution  q\\+e.  When  6  is 
zero  (Figure  III.A. 8)  it  is  seen  q  =  q\\+e  is  the  correct  root,  as  q\\  must  be  greater  than  k\\. 


frcow 


Figure  III.A.7.  LO-phonon  emission  case  2,  En-En'  >  Picolo,  “sideways”  transition  always 
allowed. 
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Figure  III. A. 8.  LO-phonon  emission  case  2,  qy  and  ky  vectors  when  6=  0. 


Equation  (III. A. 62)  becomes 
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2.  LO-Phonon  Absorption 

The  absorption  case  where  V  =  ,  ,  proceeds  in  a  similar  manner  as  for  the 

1  e  -phonon  1 

emission  case.  The  main  differences  are  Nq+\  — >  Nq  and  the  delta  function  changes  for 
the  conservation  of  energy  for  the  absorption  process.  With  these  changes,  similar  to 
equation  (III.A.47),  the  equation  for  the  electron-phonon  scattering  rate  becomes 
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Now  the  inner  function  of  the  delta  function  can  be  simplified  as  follows. 


En,+Ek,-hcDLo-En-Ek=0 


-hcow-  (X  -  En, )  +  Ek,  -Ek=  0 

-  +  Eu  +q,  ~Ek  =  0  where  hoo'a  =  fioow  -  (En,  -  En ) 


(III. A. 67) 
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Again,  0  is  the  angle  between  k  and  q  .  The  expression  may  be  rearranged  in  terms  of 
the  quadratic. 
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The  roots  to  the  quadratic  are 
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(III. A. 69) 


En—En  >  Hcolo  and  the  electron  does  not  have  sufficient  in-plane  kinetic  energy,  the 
transition  is  forbidden  (for  single  phonon  absorption).  The  subbands  are  spaced  such  that 
the  energy  difference  is  larger  than  the  LO-phonon  energy,  and  the  electron  does  not  have 
enough  in-plane  kinetic  energy,  so  the  transition  is  forbidden  (Figure  III.A. 9). 
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Figure  III.A.9.  LO-phonon  absorption  case  la,  En'-En  >  hcow  and  the  electron  does  not 
have  sufficient  in-plane  kinetic  energy,  forbidden  transition. 

LO-phonon  absorption  case  lb:  Ek  >  - hco'a  =  -[ ha>w-(En’-En )]  and  aja  <  0,  i.e., 
En'-En  >  hcow  and  the  electron  has  sufficient  in-plane  kinetic  energy,  the  transition  is 
allowed.  The  subbands  are  spaced  such  that  the  energy  difference  is  larger  than  the  LO- 
phonon  energy  and  the  electron  has  enough  in-plane  kinetic  energy,  so  the  transition  is 
allowed  assuming  no  state  blocking  (Figure  III.A.  10).  There  are  two  possible  solutions 
for  this  case,  q\\±a.  Limits  on  6  are  necessary  to  keep  the  momentum  vectors  real  valued 
(from  getting  a  negative  sign  under  the  square  root).  Also  note  that  at  resonance  the 
transition  is  allowed. 


(III.A. 73) 
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Figure  III.A.  10.  LO-phonon  absorption  case  lb,  En—En  >  ha>w  and  the  electron  has 
sufficient  in-plane  kinetic  energy,  allowed  transition. 


Equation  (III. A. 72)  becomes 
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LO-phonon  absorption  case  2:  Any  Ek  is  valid  and  <x>'a  >  0,  i.e.,  En—En  <  ha>w, 
the  transition  is  allowed.  The  subbands  are  spaced  such  that  the  energy  difference  is 
smaller  than  the  phonon  energy  so  the  transition  is  always  allowed  assuming  no  state 
blocking  (Figure  III.A.  11).  This  “sideways”  transition  has  one  solution  q\\+a  which  again 
can  be  seen  when  0  is  zero. 
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fia>Lo 


Figure  III.A.ll.  LO-phonon  absorption  case  2,  En<-En  <  hcow,  “sideways”  transition 
always  allowed. 


Equation  (III. A. 72)  becomes 
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(III. A. 75) 


LO-phonon  absorption  case  3:  Any  Ek  is  valid  and  co[,  >  0,  and  En—En  <  0,  the 
transition  is  allowed.  Absorption  from  an  upper  band  edge  to  a  lower  band  edge  is  always 
allowed  assuming  no  state  blocking  (Figure  III.A.  12).  This  “upwards”  transition  has  one 
solution  q\\+a.  Phonon  absorption  scattering  rates  (all  cases)  are  much  slower  than  phonon 
emission  scattering  rates  at  low  temperatures.  This  is  because  a  phonon  must  be  present 
to  be  absorbed.  Further  for  case  3,  absorption  via  a  large  momentum  transfer  also  reduces 
the  scattering  rate. 
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Equation  (III. A. 72)  again  takes  the  same  form  as  equation  (III. A. 75). 

3.  Mean  Scattering  Rate 

Because  there  is  a  distribution  of  carrier  energies  in  the  initial  state,  the  mean 
scattering  rate  is  calculated  by  averaging  over  the  Fermi-Dirac  distribution  of  carriers  in 
the  initial  state 


W  „ 

e  —phonon,  mean 


(III.A.76) 


where  state  blocking  has  also  been  included,  and  the  electron-phonon  scattering  rate  for 
both  emission  and  absorption  can  be  written  succinctly  as 
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(III.A.77) 


The  populations  and  electron  temperatures  can  vary  greatly  between  the  subbands,  and 
these  can  only  be  determined  by  analyzing  the  entire  QC  structure.  Modeling  the  entire 
structure  using  Monte  Carlo  simulations  will  be  discussed  later. 
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B.  ELECTRON-ELECTRON  SCATTERING 

Carrier-carrier  scattering  is  a  mechanism  that  can  be  significant,  particularly  for 
closely  spaced  energy  states.  When  the  energy  spacing  between  states  is  smaller  than  the 
LO-phonon  energy,  LO-phonon  scattering  is  forbidden  except  for  electrons  that  have 
sufficiently  high  in-plane  kinetic  energy.  Because  of  this,  the  mean  electron-phonon 
scattering  rate  will  be  reduced.  Electron-electron  scattering,  as  well  as  other  single 
electron  scattering  mechanisms  (such  as  impurity  and  interface  roughness  scattering), 
may  then  be  the  dominant  scattering  mechanisms.  Intraband  electron-electron  scattering 
along  with  LO-phonon  scattering,  must  be  modeled  in  order  to  determine  the  thermal 
electron  distributions  within  subbands.  For  these  reasons,  modeling  electron-electron 
scattering  is  important. 

The  following  diagram  (Figure  III.B.l)  shows  the  indexing  used  for  electron- 
electron  scattering,  where  the  electrons  can  be  right  or  left  handed,  and  i  and  j  are  the 
initial  states  and / and  g  are  the  final  states. 


/ 


g 

j 


Figure  III.B.l.  Feynman  diagram  for  electron-electron  scattering. 


Spin  independent  interactions  combined  with  the  Pauli  exclusion  principle  lead  to  a  tenn 
in  the  form  of  SrS2,  because  S2  =  (S1+S2)2  =  Si2+2SrS2+S22.  The  spatial  part  of  the 
wavefunction  can  be  either  symmetric  or  antisymmetric. 
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(III.B.l) 


The  total  wavefunction  is  the  product  of  the  spinor  that  describe  the  possible  spin  states 
of  the  two  electrons  and  the  spatial  part.  Since  the  electrons  are  fennions,  the  total 
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wavefunction  must  be  antisymmetric  under  simultaneous  exchange  of  both  the  spin  and 
spatial  parts.  Using  Pauli  spin  matrices 


cr  = 


"o  u 

0 

"0  -0 

,  cr3  = 

a  o8 

,  ( 7~  = 

vi  Oy 

J  °y 

vO-l 

the  spin  can  be  written  in  the  usual  way  as 
1 


S  =  — <7 
2 


(III.B.2) 


(III.B.3) 


The  permutation  operator  Pu  =  P\ 2Spin)P\ -}space)  has  eigenvalue  of  1  for  the  symmetric 
case  and  -1  for  the  antisymmetric  case.  For  a  two  electron  system,  the  eigenvalue  for  the 
SrS2  triplet  symmetric  states  is  %  and  for  the  singlet  antisymmetric  state  —A  .  This  allows 
the  permutation  operator  for  the  spin  to  be  written  as72'78 


P,  =T(1  +  4S, -S2) 


which  can  be  rearranged  to  get 
2Pn(spm)  - 1  =  4Sj  •  S2 
Since 


(III.B.4) 


(III.B.5) 


'^2  ^  ®1  °2 


(III.B.6) 


the  Pauli  matrices  are  related  to  the  exchange  operator  by  (subscripts  are  the  particle 
numbers) 

o,  o2  =2Pn(spin)  -1  (III.B.7) 

The  total  spin  and  square  of  the  total  spin  is 

S  =  i<71+i<r2  (III.B.8) 

S'S  =  f  +  i®1'®2  (III.B.9) 

From  (III.B.9)  and  (III.B.7)  the  following  is  found 

S-S  =  l  +  ^  +  ^o1  •  o2  =  P]2(spm)  +  1  (III. B. 10) 
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The  total  spin  for  the  symmetric  case  is  1  and  for  the  antisymmetric  case  0  (the 
eigenvalues  of  CT1CT2  are  1  for  the  symmetric  case  and  -3  for  the  anitsymmetric  case).  The 
total  wavefunction  may  be  written  as 

\Ys)  =  \?Ys)  and  \Va)  =  \?Va)  (IITB.ll) 

Electrons  are  indistinguishable  particles.  Denoting  the  matrix  element  as  M,  the 
following  subscript  notation  is  used. 


Mn 

=M^8f 


(III. B.  12) 


Antiparallel  spin  electrons  do  not  interfere  while  parallel  spin  electrons  do  interfere.  This 
is  because  Pauli’s  exclusion  principle  states  that  no  two  electrons  with  the  same  quantum 
numbers  (including  spin),  can  occupy  the  same  space.  Since  antiparallel  spin  electrons  do 
not  change  sign  upon  exchange,  the  square  of  matrix  element  for  the  antiparallel  spin 
case  can  be  written  as 

M2a  =\Muf +\M2lf  (III. B.  13) 

Parallel  spin  electrons  interfere  and  the  wavefunction  changes  sign  upon  exchange.  Thus, 

7Q 

the  square  of  the  matrix  element  for  the  parallel  spin  case  may  be  written  as 
M2p=\Mn  -M21 12  =  M2 -2  M12M21  (III. B.  14) 

where  the  minus  sign  comes  about  from  the  exchange.  The  matrix  element  is  the  same  as 
that  for  the  antiparallel  case,  except  for  the  additional  exchange  term  (negative  signed 
quantity).  This  exchange  effect  results  in  a  lowering  of  the  scattering  rate  for  parallel  spin 
electrons  as  compared  to  the  antiparallel  case. 

If  we  assume  equal  distributions  of  antiparallel  and  parallel  spins,  the  matrix 
element  may  first  be  written  by  including  the  Vi  probability  factors  as 

M'2  =  X-M2a  +  X-M2p  =M2-MnM21  (III. B. 15) 

Now  since  =  My-,gf,  an  additional  Vi  factor  must  be  included  to  keep  from  double 
counting. 

M 2  =  \m'2  =  |(M' -M,2U„)=  i(M,2M12  +M2IM2, -MnM2l)  (III.B  .16) 
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1.  Antiparallel  Spins 

The  first  case  considered  is  where  the  electrons  are  assumed  to  have  antiparallel 
spins,  so  no  interference  occurs  and  the  form  of  equation  (III.B.13)  applies.  All  of  the 
vectors  (k  vectors)  in  this  section  represent  parallel  component  vectors  (kj|)  and  the 
parallel  subscript  has  been  dropped  to  simplify  the  notation.  The  perturbing  potential  for 
the  electron-electron  scattering  interaction  is  the  Coulombic  potential  V  =  Vcouiomb- 
Working  in  the  Born  approximation  and  in  the  (x|  representation,  the  following  matrix 
element  is  computed  where  (jci,  X2,  X3)  =  (x,  y,  z)  with  x\,  xi  parallel  to  the  growth 
direction  X3. 


(fs\V\i,j)  =  (fg\-T^L-v\i,j) 


£X-X 


a 

~sA2 


J d 2x^d 2  x'^dx3dx3i/f(x3 )//;  (x3)i//*g(x'3)y/  j(x'3)- 


a'[(k,  -X||  +k  j  -x[|  )-(k  f  -Xu  +k g  -x[ )] 


x-x 


(III. B.  17) 


This  can  be  rewritten  as  follows,  showing  explicitly  the  X3  and  in-plane  parts  of  the 
integral  separately. 


a 

~sA2 


J dx3dx'3i//*f(x3)i//i  (x,)i//fyx()i///(x()| d2xlld2x'l 
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;[(k,-  -X,,  +k  ■  -xf,  )-(k  , -xH  +k  -xf.  ] 


Xm-Xm  +  (x-j  —  x'3 


(III.B.18) 


In  order  to  simplify  this  expression,  the  in-plane  part  which  is  essentially  a  2D  Fourier 
transform  can  be  evaluated  first.  This  can  be  done  as  follows.  The  in-plane  part  can  be 
written  as 


J  d2x^d 


ve 


'[(k,  -x,!  +k j  -xj  )-(kr  -X||  +kg  -xj )] 


Xm  -  Xn  +  (Xj  -  X3  )2 


=  J  d2x^d2X\  | 


/Ak-Xn  +/Ak'-x[, 


X|j  +  (x3  x3 ) 


(III. B.  19) 


where  Ak  =  k—  k/-  and  Ak'  =  k-kg.  Using  the  following  Bessel  function  identity1 


80,81 


the  expression  becomes 


(III.B.20) 
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J  d2X||f/2X|j 


i[(k,.-x||+ky-i[|)-(k/-x||+ki-x[|)] 


J d\d2^  e'Akx  'Ak'  x'  J J0(9||[ 
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=  J  J  J2^|J2^||eiAk-(X|l^)+iAk-xi+iAk'-xie^X3^1 
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Now  define  /?  =  X|pX|j,  and  note  that  d2x\\  =  \ch\  =  |d(x|pX|j)|  =  \dp\  since 
with  respect  to  the  |</x|||  integration,  thus  d2x\\  =  \dp\  =  d(f>pdp 

oo 

=  ||  J2X||J9||eZ(Ak+Ak)’x"e  X*^d</)pdp  e/Ak pJ0(q\\P) 
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00 

=  J  J  i/2x|,^||e'(Ak+Ak')x'  e~?" X3_x;  |  dfadpe^^JJqp) 

0 

82 

Using  the  following  equations  “ 

271 

J  d(^mpm^~"m  =  2  n  imJm(Akp ) 

0 

with  m  =  0  this  becomes 

2tc 

J  d<fciAkpcos *  =  2nJ0(Akp) 
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and 
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the  equation  now  becomes 
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(III. B. 21) 
l]  is  constant 

(III.B.22) 


(III.B.23) 

(III. B. 24) 

(III.B.25) 


(III.B.26) 
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with  q ||  =  A k  =  |k— k/|  from  evaluating  the  integral  with  the  delta  function.  Now  from  the 
d"x ||  integral  which  has  finite  limits  of  integration,  Ak  =  Ak'  is  the  only  value  such  that  the 
integral  will  be  nonvanishing,  which  yields  a  delta  function  representing  the  conservation 
of  momentum. 
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(III. B. 27) 


-9||  |*3  —  *3 


A(k/+kj;-ki.-k/) 


=  —e^l"33]AS(  Ak  +  Ak')  = 
du  / 

This  agrees  with  the  given  reference.76  This  expression  can  now  be  substituted  back  into 

the  original  matrix  element  to  arrive  at 

2  7T0cA 

I  dx-<dx\\i/f(x-<Wi(XT.)VAx\Wi(x\)Q  "v 

(III. B. 28) 


(f  s\V\i,j)  =  dx3dx'3>/f(x3)i//l  (x3V»,(*;)e 
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which  can  be  rewritten  in  the  usual  notation  as  (noting  a  =  e/4n ) 


(f.g\y\‘j)=^Am(<i,mkf+kl  k  -k,> 


(III.B.29) 


where 


Ajfg (/)  =  J dx3dx'j/f  (x3 V,.  (x3 )>// (x3 )  if/j (x3 )  e~?"  |jt3  H- 


(III. B. 30) 


and  again  q\\  =  Ak  =  | k, — k/| . 

The  scattering  rate  can  be  found  by  substituting  the  above  matrix  element  into 
Fenni’s  golden  rule  summing  over  the  final  states/ and  g.  Inserting  Pi  and  c  back  into  the 
equations 
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(III. B. 31) 

where  ET  is  total  energy,  i.e.,  the  subband  edge  plus  the  in-plane  kinetic  energy. 
Converting  the  final  states  summations  to  integrals  in  k-space 
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(III. B. 32) 
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Integrating  over  all  the  states  of  the  j  index  electron  and  including  Fermi-Dirac  functions 
to  take  into  account  the  occupancy  distributions 
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(III.B.33) 

Grouping  the  occupancy  distribution  functions  into  one  term  Fjfg,  and  evaluating  the  d2kg 
integral  to  get  rid  of  the  momentum  delta  function 
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(III.B.34) 

This  equation  gives  the  scattering  rate  of  an  electron  with  momentum  vector  k,  averaged 
over  the  other  initial  particle  states  ky.  Since  kg  can  be  written  in  terms  of  the  other  three 
momentum  vectors,  the  only  unknown  is  ky  and  this  can  be  determined  as  follows  by 

76  53 

using  relative  momentum  vectors  ’ 

k,  =k-k,  (III. B. 35) 

kfg  =  kg~kf  (III. B. 36) 
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In  this  scheme,  we  look  to  replace  the  integration  over  d2kf  with  an  integral  over  d2kfg. 
The  following  can  now  be  written  for  k / 


k/=k,+k.-kg=k,+k.-k/g-k/ 


-(k;+k.-k/g) 


(III. B. 37) 


Since  k fg  is  the  only  one  of  the  two  relative  vectors  that  is  a  function  of  k f 
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and  thus 


d"kf  =  dkf  dk r = 
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fg 


(III. B. 38) 
(III. B. 39) 

(III. B. 40) 
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as  originally  shown  in  the  given  reference.  Now  d~kfg  =  d6kfgdkfg,  where  6  is  defined  in 

53 

Figure  III.B.2,  using  notation  similar  to  the  given  reference. 

It  is  often  reasonable  to  approximate  that  the  occupancy  of  the  final  states/ and  g 
are  small,  and  ignore  state  blocking.  This  effectively  means  that  Fjfg  ~  fn\  [kj)  since  fFD / 
=  fpD,g  ~  0.  Including  state  blocking  decreases  the  scattering  rate  only  if  the  densities  in 
the  final  states  are  appreciable.  In  order  to  solve  for  q\\,  the  law  of  cosines  is  used  as 
follows 


kL  ij  =  kf  +  kj  -  2kikj  C0S<A  -  A) 
=  k 2  +  k  ■  +  Ikjkj  cos(a[) 

kL  ij  =  k}  +  K  ~  2 kfkg  cos (n  -  P) 
=  k2f  +  k2  +  2 kfkg  cos (/?) 

and 

k2  =  k2  +  k2  -  Ikjkj  cos(«, ) 

kl = kf + K  ~  2kfks  cos(A) 


(III. B. 41) 

(III. B. 42) 

(III. B. 43) 
(III. B. 44) 
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Figure  III.B.2.  Momentum  vectors  used  in  the  calculation  of  the  electron-electron 
scattering  rate. 


Summing  equations  (III. B. 41)  and  (III.B.43)  and  summing  equations  (III.B.42)  and 
(III. B. 44) 

(III. B. 45) 
(III.  B.  46) 

Substituting  equations  (III. B. 44)  and  (III. B. 45)  to  get  rid  of  k2SUmij 
k\  =  k2  +  2 (k2f  +  k2g  -  kf  -  k) )  (III. B. 47) 


k2SUmij  +  Kj  =  2k?  +  2k j 


Kumi^Ks=2k^2k\ 


The  energy  delta  function  yields  conservation  of  energy. 
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Substituting  equation  (III. B. 49)  into  (III.B.47) 


(III. B. 48) 
(III. B. 49) 
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(III.  B.  50) 
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2  2 

using  the  usual  definition  k0  =  4mltr(Ei+ErEj  -Eg).  Recall  from  the  definitions  of  cj\\,  ky, 
and  kfg,  we  may  write 

?|=|k,-k/|  =  fecfil  (III.B.51) 

Now  from  the  law  of  cosines 

(2 q^)1  =  ky  +  k2fg  -  2ky  kfg  cos  9  (III. B. 52) 

Substituting  kfg  from  equation  (III. B. 50)  into  equation  (III.B.52)  gives 

(2 q/  =  2 ky  +  k2  -Iky  -\jky  +  k2  cos0  (III.B.53) 


Performing  the  dkfg  part  of  the  d2kfg  =  d 0kf_,dkfA  integration  and  using  equation  (III. B. 47) 
gets  rid  of  the  delta  function  and  gives  q\\,  which  is  a  function  of  k/g  and  6.  Thus,  the  d6 
part  of  this  integration  remains. 
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(III.B.54) 

Now  d2kj  =  daikjdkj,  where  a\  is  the  angle  between  k,  and  k;  (the  integral  is  perfonned 
along  ky  with  angle  a\  ). 
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(III. B. 55) 


The  upper  limit  of  integration  for  kj  can  be  taken  to  be  the  top  of  the  well  barrier. 


k  =  1 2m  (v  -  F  ) 

^  j  max  T  +2  x  barrier  ^  j  ) 

V  n 


(III. B. 56) 


Equation  (III. B. 55)  is  an  expression  for  calculating  the  electron-electron  scattering  rate. 

2.  Screening 

The  formulation  so  far  has  assumed  that  there  are  only  the  two  initial  electrons 
and  they  scatter  due  to  the  Coulomb  force  between  them.  That  is,  it  has  been  assumed 
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that  there  are  no  other  electrons  (free  carriers)  present.  In  reality,  there  will  be  carrier 
concentrations  in  the  subbands.  These  free  carriers  respond  under  the  influence  of  the 
electrostatic  field.  This  will  have  the  effect  of  reducing  the  force  between  the  initial 
electrons  and  reduce  the  scattering  rate  due  to  this  screening. 

One  screening  model,  is  the  static  single  subband  screening  model,  where  the 
dielectric  function  is  modified  to  take  the  screening  into  account.  This  was  first  derived 
by  the  given  reference  and  later  used  by  others.  This  model  can  be  included  by  taking 
s  inside  the  scattering  integral,  and  replacing  it  with  a  =  sesc,  where 

^=l  +  ^-1n„(«ll)4„(«,)  (III. B. 57) 

2s  q 


and  the  polarization  factor  ignoring  collision  broadening  at  T  =  0  K  is 
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where  l%||-2%)  is  a  unit  step  function  that  when  q\\>  2 Kf  is  equal  to  1,  and  otherwise  0. 
At  any  temperature  T,  under  the  assumption  of  noninteracting  polarizability,  the 

on 

polarization  factor  becomes 
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(III. B. 59) 
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where  kn  is  Boltzman’s  constant  and  Kf  is  the  Fermi  wave  vector  defined  as 


Kf  = 


2  n  n. 


g 


(III. B. 60) 


with  the  valley  degeneracy  factor  g  =  1  for  GaAs,  which  gives  the  number  of  equivalent 
energy  bands. 

As  an  example  of  how  screening  affects  electron-electron  scattering,  consider  the 
1100  scattering  process  for  an  infinite  GaAs  well  40  nm  in  length  (Figure  III.B.3).  The 
rates  are  calculated  for  subband  populations  of  1x10  and  1x10  “  cm  at  T  =  300  K, 
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unscreened  and  including  screening.  As  expected,  the  scattering  rate  increases  as  the 
subband  population  increases.  The  screening  effect  becomes  significant  for  large  sheet 

i  o  _ o 

densities.  Sheet  densities  on  the  order  of  10  cm  “  are  two  orders  of  magnitude  larger 
than  typically  used  in  QCLs. 


Figure  III.B.3.  Electron-electron  scattering  rate  as  a  function  of  total  energy  at  T  =  300  K, 
with  and  without  including  screening,  for  subband  populations  of  (a)  lx  10 10  and  (b) 
lxlO12  cnT2. 


3.  Mean  Scattering  Rate 

The  electron-electron  scattering  rate  calculated  so  far,  is  for  the  i  indexed  electron 
at  a  given  energy  averaged  over  the  j  indexed  electron  initial  state  distribution.  Since 
there  is  a  distribution  of  carriers  in  the  initial  state,  the  mean  scattering  rate  averaged  over 
the  Fermi-Dirac  distribution  of  carriers  in  the  initial  state  can  be  found  as  follows 


[dE.W _ fFn,(k,)  m  - 

W-  _  — ‘  e~e’'  =-^-{dEkW  fFDi(ki) 

e  ~e  'mean  j dEk  fFD  i(kj )  1  e  ~e  ’ 

since 


(III. B. 61) 


(III. B. 62) 
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4.  Antiparallel  and  Parallel  Spins 

In  this  section,  taking  into  account  spin  and  state  blocking  will  be  discussed. 
Starting  with  equation  (III. B. 33),  the  solution  can  be  found  by  defining  angles  between 

on 

the  momentum  vectors  different  than  what  was  used  earlier.  Including  screening,  the 
equation  can  be  written  as 
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e  -e  32 7r3hs2 
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(III. B. 63) 


x  S(kf  +k.,  - k,.  -  k j)S(Erf  +  ETg  - ETi  -  ETj) 
rating  over  the  d2kg  to  gc 
da\kjdkj,  the  equation  becomes 


Integrating  over  the  d2kg  to  get  rid  of  the  momentum  delta  function  and  using  d2kj  = 
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From  the  usual  conservation  of  energy  and  momentum  we  have  the  following 


(III.B.64) 
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(III. B. 65) 
(III. B. 66) 


Now  defining  the  angles  for  the  momentum  vectors  as  follows  (Figure  III.B.4) 


Figure  III.B.4.  Momentum  vectors  used  in  the  calculation  of  the  electron-electron 
scattering  rate,  for  the  case  when  state  blocking  is  included. 

kg  can  be  solved  for  by  squaring  the  conservation  of  momentum  equation 
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(III. B. 67) 


K  =  (k!  +  ky  -  k/  )'  (k;  +  ky  _  k/  ) 

=  k 2  +  kj  +  k2  -  2 kjkf  cos(y )  +  2kikj  005(0;!)  -  2 kjkf  cos(y  -a,) 

Substituting  equation  (III. B. 67)  into  equation  (III.B.65) 

h2k2  fi2k2:  Trkl  f,2 

E.  + - *-  +  E .  + - J-  =  Ef+ - L  +  E+ - 

2  m  2m  2m  s  2m 

x  \k2  +  k2  +  k2  -  2kjkf  cos(y)  +  2kikj  cos(a!)  -  2 kkf  cos(y  -  a{)\ 
The  resulting  quadratic  can  be  solved  to  find  kj 


(III. B. 68) 
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kf  =  -ki  cos  (r)  +  -  k,  cos  (a,  -y)±- 
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(III. B. 69) 


From  Figure  III.B.4  the  following  equations  for  k„  and  q2  are  found,  where  angle  /is  the 
angle  between  k,-  and  k f 

kg  =  [(A;,.  +  kj  cosia^-kf  cos(/))2  +  (A:y  sin(aj) - A:/  sin(/))2]2  (III. B. 70) 

q2  =  k2  +  k)  -  2kjkf  cos  (/)  (III. B. 71) 

These  equations  are  in  a  form  so  that  state  blocking  can  be  included.  Now  to  take  into 
account  spin,  consider  the  case  where  the  electrons  have  equal  probability  of  having 
parallel  and  antiparallel  spins.  Including  the  exchange  effect  discussed  earlier,  it  is  seen 
that  the  following  substitutions  must  be  made  to  equation  (III. B. 64) 

(III. B. 72) 
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where 

q'2  =  k2  +  k2  -  2 kjkf  cos(aj  -  /) 


(III. B. 73) 


8 

and  angle  y—a\  is  between  k/and  k,.  This  agrees  with  that  in  the  given  reference. 


C.  IMPURITY  SCATTERING  AND  INTERFACE  ROUGHNESS 
SCATTERING 

In  this  section,  two  additional  scattering  mechanisms,  impurity  scattering  and 
interface  roughness  scattering,  will  be  discussed.  Beginning  with  impurity  scattering, 
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ionized  impurities  within  a  semiconductor  lattice  and  its  distributed  charges,  can  cause 
scattering  due  to  the  Coulombic  interaction.  The  corresponding  matrix  element  may  be 
written  as 


i/r  !<•)=</ 
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(III.C.l) 


This  matrix  element  can  be  simplified  for  scattering  from  an  ionized  impurity  at  X3  =  X3, 
in  a  similar  manner  as  was  done  with  the  electron-electron  scattering  case,  inserting  Pi  and 
c  back  into  the  equation  this  becomes 

(/M*)  =  dxwM?Wi  (*3)  e-9"1*3”*3 

CJ'  (III.C.2) 

e2  4(g||) 

2sAqn 

where 


AM)  =  J dx3Vf(x 3  V,.  (x3)  e^"1*3-*31  (III.C.3) 

and  q\\  =  A k  =  |k—k/(  =  k2  +kj  -2kikf  cos#  The  scattering  rate  can  be  found  by 
substituting  the  above  matrix  element  into  Fermi’s  golden  rule,  summing  over  the  final 
states 
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S(ETff-ETj) 


where  £7  is  total  energy,  i.e.,  the  subband  edge  plus  the  in-plane  kinetic  energy. 
Converting  the  final  states  summation  to  an  integral  in  k-space,  this  becomes 
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Performing  the  dkf  part  of  the  akf  =  dOkfdkf  integration  gets  rid  of  the  delta  function  and 
the  following  is  found 


AM 


impurity  %nh2s2A- 


(III.C.6) 


To  apply  this  for  the  case  of  doping  over  a  distribution,  where  a  sheet  charge  is  dx 3«3d(x3) 
at  X3  =  X3,  the  final  expression  for  computing  scattering  due  to  ionized  impurities  becomes 
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Monolayer  fluctuations  are  often  formed  at  heterostructure  interfaces.  Scattering 
can  occur  from  these  imperfect  growth  surfaces,  where  these  variations  in  the  barrier 
thicknesses  gives  way  to  variations  in  the  energy  levels  and  wavefunctions  of  the 
subbands.  Assuming  that  the  roughness  height  A(x)  at  the  in-plane  position  has  a 


correlation  function  of 


(III.C.8) 


where  the  A  is  the  roughness  and  A  is  the  correlation  length,  the  matrix  element  is  given 
by 


(III.C.9) 
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Substituting  the  matrix  element  into  Fermi’s  golden  rule  and  summing  over  the  final 
states  we  find 
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This  is  the  scattering  rate  due  to  interface  roughness  using  an  exponential  autocovariance 
function.  A  roughness  a/2,  where  a  is  the  lattice  parameter,  with  a  correlation  length  of  5 
nm  was  used  in  some  of  the  simulations  that  will  be  discussed  later. 

1.  Mean  Scattering  Rate 

Because  there  is  a  distribution  of  carrier  energies  in  the  initial  state,  the  mean 
scattering  rate  is  again  calculated  by  averaging  over  the  Fermi-Dirac  distribution  of 
carriers  in  the  initial  state. 
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D.  RADIATIVE  TRANSITION 

The  optical  transition  occurs  when  an  electron  releases  a  photon  and  transitions 
from  the  upper  to  lower  lasing  state.  As  will  be  shown,  spontaneous  emission  lifetimes 
are  on  the  order  of  microseconds,  while  other  rates  previously  discussed  such  as  phonon 
rates  can  have  lifetimes  on  the  order  of  picoseconds.  Thus,  state  lifetimes  are  not 
determined  from  spontaneous  emission  rates.  This  section  will  discuss  the  optical 
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transition  and  the  very  important  parameter  the  gain.  It  will  conclude  by  discussing  some 
fundamental  parameters  found  by  approximate  rate  equation  analysis. 

The  potential  for  the  radiative  transition  can  be  found  by  starting  with  the 
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nonrelativistic  electromagnetic  Lagrangian,  and  then  finding  the  Hamiltonian  in  the 
usual  way  H  =  ^jpiqi  -LEU . 

H  =  -1- (p2  -  2qc A  •  p  +  q) A2 )+  qcO  «  +  qcO  -  —  A  •  p  (III.D.  1 ) 

2  m  2m  m 

The  Hamiltonian  is  now  in  the  form  of  H=  H0+V.  Substituting  for  the  electron  charge,  qc 

=  —\e\,  the  interaction  is  given  by79 

V radiative  =~ AP  (III.D.2) 

m 

The  equations  for  the  E  and  B  fields  in  terms  of  the  scalar  and  vector  potentials  (O  and 
A)  can  be  written  as  a  second  rank  antisymmetric  field  strength  tensor  FMV  =  dMA  v-dvAM.H1 
These  equations  automatically  satisfy  Maxwell’s  homogeneous  equations,  and  the 
inhomogeneous  equations  are  <9 UFMV  =  JV.  The  scalar  and  vector  potentials  are  not 
uniquely  determined  because  the  following  change  AM->-Afl+dflA(x,t)  has  no  affect  on 
the  fields  E  and  B  since  8fA  v'-dvAM'  =  8MA  v-dvAM.  Therefore,  one  such  gauge 
transformation  we  are  free  to  choose  is  the  Lorentz  gauge  duAp  =  0,  which  Maxwell’s 
equations  simplify  to  one  equation  DAM  =  JjU .  Without  sources  J  =  0,  and  the  equation  for 
the  vector  potential  becomes 

V2A-^— =  0  (III.D. 3) 

dt1 

which  has  a  solution  of 


A  =  Aa  cos(q  -x-6)t)e  =  ^- +  e-i<qx-'w)]e 

Since  O  is  included  in  H0,  E  and  H  are  solved  effectively  with  0  =  0. 

5A.  zv 

E  = - =  -o)A0  sin(q  •  x  -  a>t)e 
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1  A 

H  =  —  Vx  A  = - s-sinCq-x-rufiqxe 

F  F 

The  Poynting  vector  can  now  be  found  from  these  fields 
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(III.D.7) 
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and  the  average  Poynting  vector  is  (noting  co  =  v-q) 

(S}  =  ^— q  =  ^v  (III.D.8) 

2  /j  2/u 

Inserting  Pi  and  c  back  into  the  equations,  the  average  energy  density  found  from  equation 
(II. D. 8)  is  equal  to  NPico/Y,  where  V  is  the  volume 
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from  which  Az  can  be  solved  for 
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Substituting  this  into  equation  (III.D.4),  the  vector  potential  becomes 
1  /  2PiN 
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This  can  now  be  substituted  into  equation  (III.D.2),  the  expression  for  the  radiative 
potential.  This  is  valid  for  absorption,  and  the  following  change  N  — >  N+\  must  be  made 
for  the  emission  case,  to  account  for  spontaneous  and  stimulated  emission.  In  tenns  of 
ladder  operators,  the  interaction  potential  for  the  radiative  transition  becomes  (dropping 
the  harmonic  part  and  promoting  the  momentum  term  to  the  quantum  mechanical 
operator) 
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1.  Emission  Rate 

To  find  the  2D  spontaneous  emission  rate,  equation  (III.D.12)  is  substituted  into 
Fenni’s  golden  rule  keeping  the  emission  tenn. 
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Using  the  electric  dipole  approximation,  e  iq  x  ~  1 

K.^,. = -rnr-(n + OZI*  ■  -  E,  ~ M  (hid.u) 

eV m  co  1 

Setting  N  ~  0  which  is  permissible  for  low  intensity  light,  that  is  spontaneous  emission 
occurs  without  any  photons  in  the  cavity  (to  stimulate  emission) 

~E‘  -#<o>  an  D.i5) 

£Vm  0  |  1  1 

Summing  over  the  polarizations  and  final  states  (including  all  final  photon  states  for 
spontaneous  emission  into  all  modes),  and  converting  the  k-space  summation  into  an 
integral 

(III. D.  16) 
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Note  that  d~p  =  d(j)'pdp  =  d(f)'hco/v  d(hcd) 
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To  simplify  the  bracket,  for  the  moment  we  switch  to  the  Heisenberg  picture,  where  the 
operators  (fix)  =  <f(x,t)  =  e,Ht <f(x)e~lllt  and  dx/dt  =  \/(ih)[xJ!]  =  p/m.78  Substituting  this 
commutator  into  the  bracket 

(f\p\i)  =  ij(f\[H,x\\i)  (III.D.18) 

Noting  that  the  momentum  operator  is  a  constant  of  motion,  i.e.,  p(t)  =  /;(()),  switching 
back  to  the  Schrodinger  picture  we  find 


=  (III.D.19) 

and  defining  co  =  {E-Ej)lh  (which  will  give  a  minus  sign) 

(f\p\i)  =  -imco(f\x\i)  (III.D.20) 
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Substituting  this  into  equation  (III.D.17)  and  evaluating  the  summation  for  the  2D  case, 
noting  the  polarization  sum  does  not  affect  anything  for  the  2D  case  since  the  integration 
is  only  in  terms  of  d(j)'  (i.e.,  not  d6) 
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The  final  expression  for  the  spontaneous  emission  into  all  modes  becomes  (noting 
envelope  wavefunctions  are  assumed) 
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The  3D  emission  rate  can  be  found  in  a  similar  method  as  for  the  2D  case,  except  when 
performing  the  polarization  sum  for  the  3D  case  the  square  of  the  bracket  tenn  reduces  to 
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As  can  be  seen,  there  is  a  difference  between  the  2D  and  3D  emission  cases  (noting  both 
are  single  electron  scattering  rates).  A  quick  computation  reveals  spontaneous  emission 
rates  are  on  the  order  of  microseconds,  and  thus  as  previously  mentioned,  do  not  affect 
state  lifetimes. 


2.  Gain 

The  net  transitions  from  the  upper  to  lower  lasing  state,  lead  to  an  induced  power 
=  hcoiNiWi] -M]  W\ 2)  =  hoj(ANx  W21),  where  W21  =  Wn  are  the  stimulated  emission  and 
absorption  rates  respectively.  This  amplifies  the  propagating  electromagnetic  waves 
intensity  as  dl/dx  =  gained.  The  optical  gain  of  the  medium  can  be  found,  starting  with 
equation  (III.D.14)  by  keeping  the  stimulated  emission  term.  Defining  the  photon  flux 
Iphoton  as  the  number  of  photons  per  unit  area  and  time,  substituting  for  the  momentum 
operator  and  noting  Iphotod v  =  (Number  of  photons)/V  (where  v  is  the  speed  of  light  in  the 
medium  and  V  is  the  volume),  the  stimulated  emission  rate  (a  positive  quantity)  assuming 
single  mode  becomes 
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The  gain  is  defined  as  (number  of  transitions  per  unit  volume  and  time)/(emitted  photon 
flux)  =  WILAI I  photon  =  Wl LAI piwton .  Noting  that  the  population  inversion  is  An^LA  number 
of  electrons,  the  expression  for  gain  becomes 
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Taking  into  account  a  finite  line  width,  the  delta  function  can  be  replaced  with  the 


Lorentzian 


£{hco)  = 


FWHM/2 


n  {AE  -ha>)2  +  {FWHM/2)2 


(III.D.26) 


where  FWHM  is  the  full  width  half  maximum.  The  maximum  gain  occurs  when  AE  =  hco, 
and  the  Lorentzian  is  then  equal  to  2 /{nFWHM).  Substituting  this  in  place  of  the  delta 
function,  the  expression  for  peak  gain  becomes 
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(III.D.27) 


In  Heaviside-Lorentz  units  sv  =  ninciexc  (while  in  rationalized  mksA  units  sv  =  £0nindeXc)- 
This  expression  can  be  used  for  estimating  the  threshold  population  inversion,  when  the 
threshold  gain  is  found  from  the  waveguide  resonator  analysis.  While  having  a  large 
oscillator  strength  =  2h/mco\(J\d/dx\i)\  =  2mcolh\(J\x\i)\~  is  favorable,  so  is  having  a  long 
upper  state  lifetime  to  maintain  the  3D  population  inversion  A»3d,  as  can  be  seen  from  the 
equation  for  peak  gain.  These  are  often  competing  characteristics. 


3.  Rate  Equations 

As  we  have  seen,  scattering  rates  are  density  and  temperature  dependent.  Even  to 
compute  the  mean  rates,  requires  knowing  the  densities  and  electron  temperatures  of  the 
states,  which  can  only  be  detennined  by  analyzing  the  full  structure.  Thus,  rate  equations 
are  not  well  suited  for  analyzing  full  structures.  Nevertheless,  some  important 
fundamentals  can  be  found  from  an  approximate  rate  equation  analysis.  The  rate 
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equations  for  a  three  level  system,  in  terms  of  3D  volumetric  densities,  for  the  optical 
mode  and  state  populations  are 


(n2  -nl)TYnp 


1  spon  single  mode 


T  T 

spon  all  modes  photon 


J  (ih-nJVn 
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(III.D.28) 
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where  the  upper  and  lower  lasting  states  are  indexed  by  2  and  1  respectively.  The 
notation  used  is  nphoton  is  the  3D  photon  population  mode  density,  zspon  singie  mode  is  the 
spontaneous  emission  into  the  single  lasing  mode,  rspon  aii  modes  is  the  spontaneous 
emission  into  all  modes  (noting  I  /  rspon  aU  modes  =  Nmodes  1  /%,0„  single  mode  and  the  total 
stimulated  emission  rate  1  l^stim  V^photon  ^  tspon  single  mode )?  ^photon  r(  dr,i  (So)  is  the 
photon  cavity  lifetime,  rm  =  vam  is  the  mirror  escape  time,  N  is  the  number  of  sections,  V 
=  LA  is  the  volume  of  one  section,  Vc  is  the  effective  volume  of  the  cavity,  T  =  NV/YC  is 
the  confinement  factor,  77  is  the  injection  efficiency,  and  J  =  I/A  is  the  2D  current 
density.  At  steady  state  with  z*2  ■'  tspon  all  modes  (since  Z21  '  Lpon  all  modes)  these 
equations  become 
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(III. D. 32) 
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Above  threshold,  lots  of  photons  are  present  in  the  cavity  and  Vn photon  is  large, 
and  the  second  term  on  the  right  hand  side  of  equation  (III.D.3 1)  can  be  ignored. 
(n2-n1)rV_  1  nnnui 


T  T 

spon  single  mode  photon 


(III. D. 34) 
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From  equations  (III.D.32)  and  (III. D. 33),  the  3D  population  inversion  A n3D  =  (nr-n\)  can 
be  found.  Substituting  this  back  into  equation  (III.D.34),  the  photon  density  is  determined 
to  be 


'  photon 


r  77, 


photon 


AeL 
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where  the  internal  quantum  efficiency  is 
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and  the  threshold  current  satisfying  equation  (III.D.34)  is 
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(III. D. 36) 


(III. D. 37) 
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The  power  emitted  by  the  laser  can  be  found  from  P  =  ti  (oVcnphoton /  Tm,  substituting  for 

m photon-,  t-m-i  and  r 


P  =  TJi  a,n  N^(l-Ith)  (III. D. 38) 

+  aw  N 

The  first  three  terms  in  equation  (III. D. 38)  grouped  together,  is  the  differential  quantum 
efficiency,  which  is  the  number  of  photons  emitted  per  electron.  The  internal  quantum 
efficiency  (defined  above  at  some  biases)  shows  the  output  power  is  reduced  by  an 
imperfect  (less  than  unity)  injection  efficiency  rj  and  a  finite  lower  lasing  state  lifetime 


n. 


Below  threshold,  stimulated  emission  can  be  ignored  because  there  are  few 
photons  within  the  cavity,  Vnpiw,on  ~  0.  Solving  for  n2  from  the  equation  (III.D.32)  and  for 
n\  from  the  equation  (III.D.33),  a  relation  for  the  population  inversion  is  found,  noting 
that  «2  =  rj(J!\e\L)T2.  Writing  this  in  terms  of  the  2D  population  inversion 
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(III. D. 39) 


61 


This  equation  allows  for  estimation  of  the  threshold  current.  This  simply  shows  it  is 
necessary  to  have  t\  <  r?i  for  a  population  inversion,  that  it  is  desirable  to  have  a  perfect 
injection  efficiency  77  =  I  and  a  long  upper  state  lifetime.  Improving  any  of  these 
parameters  would  serve  to  improve  the  population  inversion.  Because  rate  equations  are 
not  well  suited  for  analyzing  the  electron  transport  in  QC  structures,  in  the  next  section,  it 
will  be  discussed  how  Monte  Carlo  simulations  can  be  used  to  model  these  structures. 
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IV.  MONTE  CARLO  SIMULATION 


Since  the  lifetimes  of  the  electronic  states  within  QC  structures  are  not  constant, 
rate  equations  are  not  well  suited  for  analyzing  complete  structures  because  they  would 
be  nonlinear.  Thus,  one  approach  that  is  more  suited  for  this  type  of  transport  analysis  is 
that  based  on  Monte  Carlo  simulations.  This  approach  does  not  rely  on  any  assumptions 
about  the  electron  distributions  and  can  handle  density  and  temperature  dependent 
scattering  mechanisms.  By  analyzing  the  complete  structure,  subband  populations  and 
electron  distributions  for  all  energy  states  can  be  determined.  Quantities  such  as  current 
density  (which  is  measurable)  and  other  quantities  such  as  gain  can  be  found.  This  is 
important  because  while  having  a  large  oscillator  strength  is  favorable,  so  is  having  a 
long  upper  state  lifetime  to  maintain  the  3D  population  inversion  A/?3D,  as  can  be  seen 
from  the  equation  for  peak  gain  =  2e  h  An  D/  m  ~  s  v  coF  W 1 1 M\ (J\dl dx  |  /)  |  (equation 
(III.D.27)).  These  are  often  competing  characteristics. 

A.  ELECTRON  TRANSPORT  IN  QC  STRUCTURES 

As  has  been  discussed,  the  electron  transport  in  entire  QC  structures  is  not  suited 
for  rate  equation  analysis,  and  the  large  number  of  scattering  rates  that  must  be  computed 
would  further  complicate  attempting  to  use  such  an  approach.  From  the  transport  analysis 
of  QC  structures,  we  ultimately  seek  to  find  the  populations  and  electron  distributions  of 
all  states  in  the  structure.  The  time  evolution  of  the  electron  distribution  functions  is 
described  by 

%  =  ElvwU 1  -  /J-  ■W,'/w( i  -  fvi  (iv.a.1) 

at  k>p> 

where  S  is  the  total  scattering  rate  (sum  of  all  mechanisms).  The  Monte  Carlo  method, 
applied  to  analyzing  the  electron  transport  in  QC  structures,  consists  of  tracing  a 
sufficient  number  of  electrons  throughout  the  structure  over  time,  such  that  the  results 
will  be  a  good  representation  of  the  electron  ensemble  within  the  device.  It  relies  on  the 
use  of  random  number  generation,  to  appropriately  choose  from  the  probability 
distributions  encountered  in  this  type  of  a  simulation.  There  are  several  equivalent  ways 
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to  formulate  this  technique.95'99  In  what  follows,  it  will  be  discussed  how  we  chose  to 
implement  this  technique. 


B.  MONTE  CARLO  METHOD 

In  general,  the  motion  of  a  charged  particle  is  due  to  the  presence  of  fields, 
electric  and  magnetic,  and  scattering  events.  The  time  durations  between  scattering 
events,  the  scattering  times,  historically  has  been  called  the  free  flight  times.  During  its 
time  evolution,  the  electrons  subband  and  momentum  change,  and  so  does  the  total 
scattering  rate  for  all  mechanisms.  Consequently,  the  free  flight  times  will  change  too.  In 
order  to  be  able  to  trace  electrons  throughout  a  structure,  we  must  develop  a  method  for 
stochastically  generating  these  flight  times.  For  example,  consider  the  simple  case  of  an 
electron  moving  under  the  influence  of  an  applied  electric  field  in  the  Xi  direction,  the  k\ 
component  changes  as  ki=ki>t=o-\e\/h\Ei\t,  while  the  other  momentum  vector  components 
remain  unchanged.  At  the  end  of  a  free  flight,  the  electrons  momentum  and  energy  are 
updated,  and  the  electron  is  then  scattering  into  its  next  state. 

In  QC  structures,  an  electron  is  characterized  by  its  subband  and  components  of 
the  in-plane  momentum  vector  k.  Because  the  wavefunctions  already  contain  the  effects 
of  the  applied  electric  field,  there  is  no  need  to  include  the  momentum  component  that  is 
in  the  growth  direction.  Next  it  will  be  discussed  how  to  stochastically  generate  these  free 
flight  times  for  the  special  case  when  the  scattering  rate  is  constant.  As  will  be  shown, 
this  fonnulation  is  still  applicable  to  the  more  practical  case  where  the  scattering  rate 
changes  over  the  simulation  time,  as  will  be  necessary  for  modeling  the  electrons  within  a 
QC  structure. 

1.  Flight  Times 

Considering  the  scattering  rate  out  of  a  subband  to  be  constant  ro=  1  It,  an 
electrons  probability  to  undergo  a  scattering  event  in  a  duration  dt  is  T0dt,  and  dn  = 
-T0ndt.  This  simple  equation  can  be  solved  to  find  how  the  population  changes  with  time. 
n(t)  =  n(  0)e“r»'  (IV.B.l) 

The  probability  that  the  electron  will  scatter  between  t  and  t+dt,  depends  on  the  number 
of  electrons  in  the  subband  and  the  scattering  rate.  This  relation  may  be  written  as 
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dtP(t)  =  dtT0e~rj  (IV.B.2) 

where  P  is  a  random  number  from  0  to  1.  Using  the  direct  continuous  technique97  to 
stochastically  generate  the  flight  times,  from  equation  (IV.B.2) 

'/ 

r  =  |  dtP(t)  (IV.B.3) 

0 

where  r  is  a  random  number  between  0  and  1 .  The  solution  for  this  equation  is 
r  =  -Q~v°tf  + 1  (IV.B.4) 

from  which  the  flight  time  tf  =  -1/  roln(l-r)  is  solve  for.  Since  r  is  a  number  between  0 
and  1,  so  is  (1  -r),  and  it  is  equivalent  to  write 

tf=-^rHr)  (IV.B.5) 

o 

This  expression  allows  us  to  stochastically  generate  the  flight  times  for  a  constant  total 
scattering  rate. 

As  discussed  earlier,  the  scattering  rates  will  not  be  constant.  Nevertheless,  this 
formulation  is  still  applicable  to  the  case  where  the  total  scattering  rate  changes  over  the 
simulation  time.  To  deal  with  this,  we  will  use  a  clever  technique  formulated  in  the  given 
references,95'96  called  self  scattering.  Self  scattering  is  the  difference  between  the  constant 
scattering  rate  and  the  actual  rate.  The  actual  scattering  rate  is  not  constant  and  will 
depend  on,  for  instance,  the  scattering  mechanism,  in-plane  momentum  vector  k,  subband 
of  the  initial  electronic  state,  electron  densities,  and  temperatures. 

ro  =Tseif  +T(k)  (IV. B. 6) 

The  reason  for  the  inclusion  of  this  self  term,  is  that  it  simplifies  the  calculations  that 
follow  by  allowing  us  to  treat  the  scattering  rate  as  being  constant.  The  self  term  will 
eventually  be  removed,  so  that  it  will  have  no  affect  on  the  final  results  of  the  simulation. 
For  this  technique  to  work,  T0  must  always  be  chosen  greater  than  the  maximum  actual 
scattering  rate  T(k).  The  penalty  paid,  is  that  at  the  end  of  every  free  flight  time,  it  must 
be  checked  to  see  if  an  actual  scattering  event  or  a  fictitious  self  scattering  event 
occurred.  If  a  fictitious  self  scattering  event  occurred,  no  real  scattering  event  occurred 
and  the  electrons  momentum  is  updated  and  a  new  flight  time  is  generated.  If  on  the  other 
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hand  it  is  detennined  that  a  real  scattering  event  occurred,  the  electrons  energy  and 
momentum  are  updated  and  it  is  then  scattered  into  its  new  state,  where  its  energy  and 
momentum  are  updated  for  this  new  state.  This  process  of  checking  to  see  if  a  real 
scattering  or  self  scattering  event  occurred,  can  generate  a  considerable  amount  of 
oversampling  (self  scattering  events).  Despite  the  oversampling  that  occurs  due  to  the 
inclusion  of  the  self  tenn,  it  is  much  better  than  the  alternative  which  would  require 

go 

changing  equation  (IV.B.2)  to  following  integral  equation  (which  is  in  the  exponent). 

t 

-JV'r(k,/') 

dtP(t)  =  dtT(k,t)e  0  (IV.B.7) 

This  equation  is  not  suited  for  efficient  repeated  numerical  calculations,  and  for  this 
reason  the  self  scattering  technique  is  used  instead. 

To  use  equation  (IV. B. 5)  to  generate  the  flight  times,  along  with  the  self 
scattering  technique,  T0  must  be  found  for  the  simulation.  There  are  different  ways  to 
implement  this,  we  choose  to  use  the  maximum  rates  defined  below,  summing  over  all 
mechanisms  (denoted  by  m)  the  maximum  scattering  rates  TmaXj  m 

r„  =Zrm..,»  (IV.B.8) 

m 

where 

rmax,m  =  max,  max, k  [r  ^y  (k)]  (IV.B.9) 

/ 

for  our  initial  mapping  scheme  in  the  selection  process  that  will  follow.  This  will  thereby 
introduce  an  additional  self  tenn  in  the  selection  process,100  that  will  also  have  to  be 
removed  using  a  straight  forward  rejection  technique.  This  provides  a  convenient  way  of 
implementing  simulations  with  a  large  number  of  scattering  rates,  with  the  alternative 
method  being  to  choose  the  maximum  rate  as  a  function  of  k  value.  This  method  of 
finding  T0  is  illustrated  graphically  in  Figure  IV.B.l. 
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Figure  IV.B.l.  Graphical  illustration  for  finding  ro,  for  the  Monte  Carlo  method. 

2.  Selection  Process 

Now  that  r0  can  be  found  and  the  flight  times  can  be  generated,  we  must  proceed 
with  the  selection  process  that  determines  the  scattering  mechanism  (or  if  a  self  scattering 
event  occurred),  the  final  state,  and  final  energy  and  momentum  of  the  followed  electron. 
Single  electron  scattering  mechanisms  are  all  handled  very  similarly,  and  electron- 
electron  scattering  is  a  special  case  requiring  more  parameters  to  be  determined.  It  will 
first  be  describe  how  the  selection  process  works  for  the  single  electron  scattering  cases, 
and  then  for  the  electron-electron  scattering  case.  Beginning  with  a  new  free  flight,  a 
random  number  r0  is  generated  between  0  and  T0  and  mapped  onto  the  different 
maximum  scattering  mechanisms  rates,  using  the  following  bin  walling  procedure,97'101 
where  the  mechanism  m  =  x  is  selected  if 

Zr— «.1  Sro  <Zr-*MM  (IV.B.10) 

m= 0  m= 0 

If  no  self  scattering  occurs  and  mechanism  m  =  x  is  selected,  the  final  subband  stat e/=  v 
is  chosen  by  generating  another  random  number  r\  between  0  and  rmax  Jn=xj  and  mapping 
as  follows 

X Fmax, *,,->/  ^  h  <2Xx, m.^f  (IV -B • 1  1 ) 

/= 0  /= 0 
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Once  the  mechanism  and  final  state  subband  have  been  determined,  the  final  momentum 
vector  kf  =  kj{ cos6*xi+  sin^x?)  must  be  found.  As  an  example,  the  magnitude  kf  is  easily 
found  from  the  conservation  of  energy  Etj  =  Et/  or  for  LO-phonon  scattering  Ej,  = 
Etj±ELo,  where  ET  =  E„+h  k  12m.  The  angle  6  can  be  found  by  using  another  Monte 
Carlo  technique.  First  the  angle  is  randomly  chosen  from  a  uniform  distribution  between 
0  and  2n.  The  self  tenn  is  removed  using  a  rejection  technique,  which  is  an  iterative 
scheme  keeping  the  choice  if 

r<^~  (IV.B.12) 

P 

max 


and  otherwise  rejecting  it.  P{9)  is  found  from  the  corresponding  scattering  mechanism 
equations  by 


P(0)  = 


K 

\dOK 


(IV.B.13) 


where  K  is  the  kernel  inside  the  scattering  rate  integral.  This  selection  process  and  how  it 
unfolds  is  graphically  illustrated  in  Figure  IV.B.2. 
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Figure  IV.B.2.  Graphical  illustration  of  the  Monte  Carlo  selection  process. 

For  electron-electron  scattering,  many  more  additional  terms  must  be  found.  The 
subband  i,  k,  and  6,  of  the  initial  followed  electron  (three  parameters)  are  known,  and 
these  three  parameters  for  the  second  or  parent  electron  to  scatter  with  must  be  found,  as 
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well  as  for  the  two  final  state  electrons.  In  our  implementation,  this  is  done  by  bin  walling 
in  a  similar  fashion  as  discussed  earlier,  to  find  the  ijfg  process  (a  mapping  look  up 
scheme  is  used).  A  parent  electron  is  then  randomly  chosen  from  the  other  electrons 
within  subband  j  of  the  ensemble,  which  determines  j,  kj,  and  Oj.  Since / and  g  are  known, 
only  four  parameters  remain  to  be  found  (kj.  Of,  kg,  and  0g).  Only  one  of  these  really  needs 
to  be  found  though,  as  the  others  are  dependant  by  the  conservation  of  energy  and 
momentum.  We  chose  ^randomly  from  a  uniform  distribution  between  0  and  2n,  then 
using  equations  from  section  III.B  (equations  (III.B.69)  et  cetera),  the  other  three 
parameters  can  be  found.  The  self  term  can  be  removed  by  solving  for  q\\  from  equations 
(III.B. 52)  or  (III.B. 71),  and  using  the  rejection  technique  of  equation  (IV.B.13). 

C.  IMPLEMENTATION 

A  flow  chart  describing  the  Monte  Carlo  simulation  is  shown  in  Figure  IV.C.l. 
The  algorithm  consists  of  calculating  the  scattering  rates  and  the  main  Monte  Carlo 
selection  process  previously  described.  Because  a  QC  structure  is  periodic  and  charge  is 
conserved,  when  an  electron  with  in-plane  wavevector  k  is  scattered  into  another  section, 
an  electron  is  introduced  with  the  same  k  vector  into  an  equivalent  state  in  the  section 
being  modeled.  Typically,  the  center  section  of  a  three  section  structure  is  modeled. 
Initially  the  electrons  are  distributed  thermally  at  the  lattice  temperature  according  to 
Fenni-Dirac  statistics,  in  one  or  more  energy  levels.  In  theory,  the  scattering  rates  need  to 
be  recalculated  every  time  an  electron  scatters  to  a  new  state.  However,  in  practical 
implementation  of  this  algorithm,  the  scattering  rates  need  only  be  calculated  after  an 
evolved  time  step  duration  At,  and  this  duration  of  time  is  chosen  such  that  the  scattering 
rates  do  not  change  too  much.  For  example,  a  At  ~  0.1  psec  might  be  needed  for  quickly 
changed  transient  times.  So  in  practice,  each  electron  in  the  ensemble  is  traced  for  a  time 
step  duration  At,  then  the  scattering  rates  are  recalculated  based  on  the  new  subband 
populations  and  electron  distributions.  In  this  way,  the  scattering  rates  are  consistently 
solved  for  using  the  correct  subband  distributions,  as  the  time  evolution  of  the  ensemble 
unfolds.  This  allows  for  state  blocking  to  be  taken  into  account,  and  ensures  that  density 
and  temperature  dependant  scattering  rates  are  correctly  modeled.  The  number  of 
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particles  must  be  chosen  enough  for  accurate  calculation,  and  typically  ~  104  is  adequate 
for  most  simulations.  Note  that  it  is  not  that  important  how  many  electrons  or  the  sheet 
density  each  particle  represents,  only  that  there  are  enough  for  accurately  statistically 
modeling  the  ensemble  within  the  structure. 


Ensemble 

loop 


Figure  IV.C.l.  A  flow  chart  for  the  Monte  Carlo  simulation  of  QC  structures. 
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Because  the  scattering  rate  calculations  are  the  time  consuming  part  of  the 
simulation,  we  choose  to  calculate  and  store  the  rates  before  each  time  duration.  That 
way,  the  rates  are  only  calculated  once  for  each  time  step  duration.  The  Monte  Carlo 
algorithm  is  continued  until  it  is  detennined  that  state  populations  have  reached  steady 
state  and  are  not  changing.  This  is  done  by  making  sure  the  variance  of  the  populations  is 
below  some  minimum  acceptable  value.  Once  steady  state  is  reached,  the  electron 
densities  and  temperatures  are  output,  and  other  parameters  such  as  various  scattering 
rate  lifetimes,  current  density,  and  gain  are  easily  found  (all  scattering  events  are  stored). 
As  an  example  of  a  typical  simulation,  shown  in  Figure  IV.C.2  is  a  Monte  Carlo 
simulation  of  QC  structure  D619F10E  which  has  six  energy  levels  within  one  section  (the 
details  of  this  structure  will  be  discussed  later),  where  all  electrons  were  initially 
distributed  in  state  2  at  the  lattice  temperature  of  T=  25  K. 


Figure  IV.C.2.  Monte  Carlo  simulation  of  structure  D619F10E  at  a  bias  of  53.6 
mV/section,  with  all  electrons  initially  distributed  in  state  2.  Inset:  Two  sections  of  the 
structure  are  shown,  with  one  section  outlined. 

Under  the  transport  picture  described  by  this  type  of  Monte  Carlo  simulation, 
where  the  scattering  between  electronic  states  is  determined  by  the  spatial  wavefunctions 
solved  via  Schrodinger’s  equation,  coherence  effects  are  not  taken  into  account.  For 
example  the  peak  current  through  a  thickened  barrier  is  not  affected,  only  its  transmission 
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sharpness  decreases,  and  this  is  not  accurate.  This  implies  that  the  scattering  is  not  a  fully 
coherent  process  and  scattering  between  weakly  coupled  states  with  a  small  anticrossing 
is  mostly  a  noncoherent  process,  due  to  dephasing  scattering  that  effectively  interrupts  or 
scrambles  the  phase  coherence.  No  phenomenological  dephasing  parameter  is  introduced 
in  our  simulations.  This  can  sometimes  lead  to  optimistic  results  as  well  as  an 
overestimation  of  the  scattering  between  states  with  a  small  anticrossing,  and  hence  an 
overestimation  of  the  current  density  when  aligned  with  weakly  coupled  parasitic  current 
channels.  The  scattering  between  more  than  one  section  is  taken  into  account  in  step  well 
structures  because  of  the  nature  of  the  step  well  QC  structure  (which  potentially  can 
further  lead  to  an  overestimation  of  the  scattering  between  weakly  coupled  states).  The 
scattering  mechanisms  that  are  included  are  LO-phonon,  electron-electron,  impurity,  and 
interface  roughness  scattering.  These  scattering  rates  are  calculated  via  Fermi’s  golden 
rule.  All  rates  from  these  simulations  are  net  rates  and  include  backscattering.  State 
blocking  and  screening  are  also  included  in  the  simulations.  Next  we  discuss  waveguide 
structures  commonly  used  in  QCLs  for  the  confinement  of  emitted  photons. 
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V.  WAVEGUIDE  AND  RESONATOR 


QCL  devices  are  fabricated  by  growing  N  repeated  sections  for  the  active  region, 
and  then  processing  the  wafer  sample  into  a  waveguide  to  form  the  laser  resonator.  In 
these  sections,  waveguide  parameters  of  interest  and  the  two  common  types  of  waveguide 
structures,  surface  plasmon  and  metal-metal  waveguides,  are  discussed.  A  third  structure, 
a  metal-patterned  substrate  waveguide  was  also  investigated. 


A.  QCL  RESONATORS 

The  waveguide  and  resonator  are  important,  because  the  losses  essentially  set  the 
threshold  condition.  At  threshold,  the  gain  equals  the  losses.  The  threshold  gain  can  be 
found  by  considering  a  one  round  trip  inside  the  resonator  as  follows 


-2  a„,L 


/  e 


gmnth 


RqR\q 


2FxgainxL  _ 


1 


1 


ln(f?„f?i)  +  aH 

V  2L 


(a  +  a  ) 

_  V  m  u'  / 


(V.1) 


where  I0  is  the  initial  intensity,  am  =  -l/2Zln(RoRi)  is  the  mirror  loss,  aw  is  the  waveguide 
loss,  L  is  the  length  of  the  resonator,  and  R  is  the  intensity  reflectivity  at  a  facet.  The 
approximation  that  the  mode  gain  =  Txgain  has  been  used.  The  confinement  factor  T,  is 
defined  as  the  ratio  of  the  guided  energy  inside  the  active  region  to  the  total  energy  of  the 
wave 
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where  S  =  'AExH*. 

Since  the  waveguide  mode  is  attenuated,  the  propagation  constant  is  complex 
/?=  j3'+ij3",  noting  that  aw  =  2(3".  In  terms  of  the  electric  field  and  considering  the  case 
where  the  field  facet  reflectivity  r  =  is  the  same  for  both  facets,  the  following  can 
be  written 
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In  general,  it  is  seen  from  the  phase  terms  that  the  Fabry-Perot  resonance  condition  is 
found  P'L  =  nn-(j).  By  ignoring  (f>,  which  is  often  done  for  optical  systems,  at  resonance 
the  usual  condition  is  found  L  =  haJ2.  Next,  the  types  of  waveguides  used  for  QCLs  will 
be  discussed. 


B.  SURFACE  PLASMON  AND  METAL-METAL  WAVEGUIDES 

Two  types  of  optical  confinement  structures  for  THz  QCLs  have  been  used,  the 
semi-insulating  (SI)  surface  plasmon  and  the  metal-metal  waveguides.1,32’102'104  A  third 
type  of  waveguide  similar  to  the  metal-metal  waveguide,  the  metal-patterned  substrate 
waveguide,  was  also  analyzed  (Figure  V.B.l). 


(a)  (b)  (c) 


Figure  V.B.l.  QCL  waveguides,  the  (a)  surface  plasmon,  (b)  metal-metal,  and  (c)  metal- 
patterned  substrate  waveguides. 

The  surface  plasmon  waveguide  (Figure  V.B.l(a))  consists  of  a  metallic  top  and 
two  side  contacts,  and  a  thin  n+  doped  semiconductor  layer  or  plasma  layer  (where  the 
real  part  of  the  permittivity  is  negative)  sandwiched  between  the  active  region  and  the  SI 
GaAs  substrate.  Gaps  are  present  between  the  sides  of  the  active  region  and  the  metallic 
lateral  contacts.  This  type  of  structure  is  easier  to  fabricate  than  the  metal-metal 
waveguide  since  it  does  not  require  substrate  removal  and  wafer  bonding.  Since  the  fields 
extend  substantially  below  and  outside  the  active  region,  the  confinement  factor  is  less 
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than  that  of  the  metal-metal  waveguides.  However,  the  advantage  can  be  improved  output 

1 3 

coupling  for  higher  power  applications. 

The  metal-metal  waveguide  consists  of  two  metallic  layers  above  and  below  the 
active  region  (Figure  V.B.l(b))  which  is  above  a  n+  GaAs  substrate.  The  fields  of  the 
mode  are  usually  well  confined  with  often  a  near  unity  confinement  factor,  which  can 
lead  to  lower  threshold  gains.  The  waveguide  width  may  often  be  made  smaller,  than  the 
surface  plasmon  waveguide  counterpart,  with  acceptable  performance.  This  can  allow  for 
reduction  in  the  device  area.  The  metal-metal  waveguide  is  also  better  suited  for  higher 
operating  temperatures.  The  disadvantage  is  that  the  output  coupling  is  worse  due  to  the 
high  reflectivity  of  the  subwavelength  rectangular  aperture,  and  the  radiated  beam  pattern 
is  much  more  fan  shaped.  At  THz  frequencies,  we  are  in  a  waveguide  region  in  which  the 
Fresnel  reflection  coefficient  is  not  valid  to  use  for  the  reflectivity  for  the  metal-metal 
waveguide.  The  high  reflectivity  is  an  advantage  in  the  sense  that  a  high  reflectivity  (HR) 
coating  on  one  facet  end  is  rarely  needed  with  this  type  of  waveguide.  Although  perhaps 
difficult  to  fabricate,  others  have  shown  improved  output  power  in  metal-metal 
waveguides  by  micromachining  a  matching  horn  antenna  on  the  facet  end.105  Both  the 
surface  plasmon  and  metal-metal  configurations  use  surface  plasmons  attached  to  the 
contact  layers. 

The  waveguide  losses  and  reflectivities  for  the  waveguides,  are  found  using  finite 
element  method  (FEM)  solvers.  '  The  Drude  model  is  used  to  determine  the 
material  parameters,  in  which  the  complex  permittivity  can  be  found  from 
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where  \e\  is  the  charge  of  the  electron,  m  is  the  effective  mass  of  the  electron,  and  co  is  the 
frequency.  In  our  simulation  model,  the  active  region  was  nominally  taken  as  GaAs  and 
Sb  was  set  to  12.4  corresponding  to  the  f=0K  value.  The  relaxation  times  r,  were  set  to 
0.1  ps  for  the  n+  layers,  0.5  ps  for  the  active  region,  and  0.06  ps  for  the  metallic  gold 
layers.108’110  For  the  active  region  a  background  carrier  concentration  of  2xl01:i  cnT3  was 
assumed,  and  a  concentration  of  5x10“'  cm  was  used  for  the  metallic  gold  layers.  The 
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plasma  and  contact  layers  concentrations  vary  by  design,  and  are  listed  in  those  particular 
sections.  From  Ampere’s  law  V  *H  =  J+dD/dt  =  ctE— 10D  and  equation  (V.4),  it  is  seen 
the  complex  permittivity  can  alternatively  be  entered  as  a  real  valued  permittivity  Sb  and  a 
complex  conductivity  a  =  n^e  r/m(l-icoT).  Additionally,  it  can  also  be  entered  as  a  real 
valued  effective  permittivity  s'  and  a  real  valued  effective  conductivity  ox"  since  V  xH  = 

ox  "E  -ico sE,  noting  that  s'  can  be  a  positive  or  negative  value.  By  considering  co  in  the 

2  1/2 

high  frequency  limit,  the  plasma  frequency  may  be  defined  as  cop  =  {nyae~lz\,m  )  .  Since 
ck=  (o)  -cop)  it  is  seen  that  when  co<  cop,  k  is  purely  imaginary  and  waves  incident  on 
the  plasma  are  reflected  and  waves  inside  fall  of  exponentially  with  distance,  thus  the 
medium  behaves  more  like  a  “metal.”  For  o)>  0)p  transmission  occurs  and  the  medium 
has  more  of  a  lossy  “dielectric”  characteristic. 

A  waveguide  similar  to  the  metal-metal  waveguide,  the  metal-patterned  substrate 
waveguide,  was  also  briefly  investigated  (Figure  V.B.l(c)).  The  idea  was  if  MBE  growth 
of  the  active  region  could  be  accomplished  on  top  of  a  patterned  substrate,  this  could 
avoid  having  to  wafer  bond  while  at  the  same  time  keep  a  structure  capable  of 
maintaining  a  high  confinement  factor.  It  was  ultimately  determined  that  the  defects 
introduced  by  the  growth  process  over  the  patterned  substrate,  would  likely  be 
unacceptable  for  high  quality  MBE  growth.  Nevertheless,  the  results  of  the  analysis  will 
be  briefly  mentioned,  as  similar  type  of  3D  analysis  could  be  useful  for  mode  suppression 
waveguides.  In  such  a  waveguide,  a  metallic  layer  would  be  deposited  on  the  top  surface 
of  the  active  region,  and  the  pattern  would  also  have  to  be  dense  enough  to  keep  the 
modes  confined  while  also  not  substantially  increasing  the  waveguide  loss.  Our  analysis 
showed  for  the  guides  we  modeled,  compared  to  metal-metal  waveguides,  that  the 
confinement  factor  was  affected  little  by  square  patterns  from  1  by  5  fjm  (1  fjm  thick 
strips  with  5  /an  spacing  between  the  strips)  to  2  by  2  /an,  but  the  waveguide  loss  was 
increased  for  the  less  dense  patterns.  It  was  also  shown  for  a  2  by  2  fjm  pattern 
waveguide,  that  the  waveguide  loss  could  be  made  lower  using  a  n+  substrate  (as 
opposed  to  a  SI  substrate  with  a  thin  n+  layer)  for  frequencies  beyond  about  2.5  to  5  THz. 
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VI.  THz  LO-PHONON  QC  STRUCTURES 


Two  different  LO-phonon  approaches  were  investigated,  the  new  step  well  QC 
structures  we  proposed,40’111'113  and  a  conventional  square  well  design.  Approximate  rate 
equation  analysis  was  initially  performed  on  the  step  well  structures,  and  the  structures 
were  analyzed  using  Monte  Carlo  simulations.  The  analysis  indicated  the  injector  of  one 
of  the  step  well  structures  could  further  be  improved,  and  simulations  on  a  modified 
structure  showed  a  more  uniform  gain  over  a  range  of  biases.  Simulations  were  also 
performed  on  a  high  power  square  well  LO-phonon  design.  Our  analysis  was  compared 
to  experimental  measurements,  and  found  to  be  in  reasonable  agreement. 

A.  STEP  WELL  QC  STRUCTURES 

Even  though  it  is  possible  to  arrange  the  energy  levels  in  a  step  well  such  that  the 
radiative  transition  and  LO-phonon  transition  are  within  the  same  step  well,  this  does  not 
mean  it  will  necessarily  be  useful  for  a  laser  application.  For  the  device  to  lase,  a 
population  inversion  must  be  kept  for  sufficient  gain  to  overcome  the  losses.  While 
intrawell  radiative  transitions  can  have  large  overlap  of  the  electron  wavefunctions  that 
yield  large  oscillator  strengths,  there  can  be  a  trade-off  between  the  oscillator  strength 
and  upper  state  lifetime,  as  the  scattering  from  the  upper  radiative  state  to  the  lower  state 
is  usually  increased  with  increasing  oscillator  strength  in  three-level  LO-phonon 
structures.  This  increase  in  oscillator  strength  is  one  of  the  reasons  for  studying  these  step 
well  structures,  as  increasing  the  oscillator  strength  can  increase  the  gain  and  hence  the 
output  power  of  a  device.  As  will  be  shown,  the  section  length  will  remain  about  the 
same  as  other  LO-phonon  designs. 

Initially,  an  approximate  rate  equation  analysis  was  perfonned  on  the  step  well 
structure  shown  previously  in  Figure  II.B.  1(b)  40111  The  scattering  rates  (W2\  and  W \ o) 
must  be  calculated  for  the  step  quantum  well  structures  in  order  to  estimate  whether  a 
population  inversion  is  likely  between  states  2  and  1 .  The  simplest  rate  equation  analysis 
allows  us  to  write  the  population  of  the  middle  state  1  (n\  )  as 
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-L  =  n2W2l-nlWl0  (VI.A.1) 

at 

where  H2  is  the  population  of  state  2.  The  necessary  and  intuitive  condition  at  steady  state 
for  a  population  inversion  to  exist  ( n2  >  n\)  is  W\o  >  W2\  (in  terms  of  corresponding 
lifetimes  rio  <  t2\).  In  general,  the  scattering  rates  are  a  combination  of  all  possible 
scattering  mechanisms,  i.e.,  electron-phonon,  electron-electron,  impurity,  and  interface 
roughness  scattering.  In  this  approximate  rate  equation  analysis,  we  estimate  the 
scattering  rates  by  taking  into  account  the  LO-phonon  and  electron-electron  scattering 
rates.  As  will  be  seen  later  from  our  Monte  Carlo  analysis,  the  electron  temperatures  can 
be  ~  50  to  100  K  or  higher  than  the  lattice  temperature.  For  the  rate  calculations  that 
follow,  the  lattice  temperature  will  be  taken  to  be  25  K  and  the  electron  temperature  to  be 
100  K. 

Figure  VI.A.l  shows  the  mean  LO-phonon  scattering  rates  for  the  4.3  THz  step 
quantum  well  previously  shown  in  Figure  II.B.l(b)  (and  shown  again  inset),  calculated 
for  initial  state  populations  from  1*109  to  lxlO10  cm-2  for  Tiattice  =  25  K  and  Te/ectmn  =  100 
K.  All  of  the  scattering  rates  are  essentially  constant  and  independent  across  the 
populations  computed.  As  expected,  the  2  to  1  transition  has  the  smallest  mean  scattering 
rate  and  hence  the  longest  lifetime,  because  the  energy  spacing  of  17.9  meV  between 
states  2  and  1  is  below  the  LO-phonon  energy  spacing  of  about  36  meV.  Only  the 
electrons  with  sufficiently  high  in-plane  kinetic  energy  (~  36-17.9  meV  =  18.1  meV  or 
higher)  can  participate  in  the  scattering  process,  while  others  are  forbidden.  The 
maximum  scattering  rate,  which  occurs  when  the  electron  has  just  enough  in-plane 
energy  such  that  the  energy  difference  is  exactly  at  the  resonant  LO-phonon  energy  Hcolo, 
is  also  shown  in  Figure  VI.A.l.  It  can  be  seen  that  the  mean  scattering  rates  for  the  1  to  0 
transition  are  faster  than  the  2  to  1  transition,  and  are  even  faster  than  the  2  to  1 
maximum  rate.  Also  included  is  the  parasitic  mean  scattering  rate  for  the  2  to  0  transition. 
This  parasitic  lifetime  is  1.8  picoseconds  (psec),  which  is  relatively  short.  These 
scattering  rates  show  that  the  mean  lifetimes  are  such  that  t2\  =  5.6  psec  >  rio  =  0.37 
psec.  This  indicates  that  the  LO-phonon  rates  are  favorable  for  keeping  a  population 
inversion.  Of  course,  the  population  inversion  is  affected  by  all  scattering  mechanisms, 
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not  just  by  the  electron-phonon  scattering  rates.  Next  we  will  calculate  the  electron- 
electron  scattering  rates  and  see  how  incorporating  these  rates  may  affect  the  population 
inversion  analysis. 


2  4  6  8  10 

n(109  cm'2) 


Figure  VI.A.l.  Mean  LO-phonon  scattering  rates  for  the  4.3  TFIz  step  quantum  well 
(shown  inset)  as  a  function  of  initial  state  populations  for  Tiattice  =  25  K  and  Teiectron  =  100 
K  (solid  lines).  For  the  2  to  1  transition  the  maximum  scattering  rate  (which  is  not  an 
averaged  rate)  is  also  shown  (dashed  line).  The  lifetimes  are  listed  (reciprocals  of  the 
essentially  constant  scattering  rates). 

The  electron-electron  scattering  rates  were  calculated  for  the  same  4.3  THz 
structure  as  for  the  LO-phonon  analysis  using  Teiectron  =  100  K  for  initial  state  populations 
of  lxlO9  to  lxlO10  cm2.  Figure  VI.A.2  shows  the  electron-electron  scattering  rates  for 
the  2  to  1  and  1  to  0  transitions  that  pertain  to  the  simplified  rate  equation  analysis 
already  discussed.  It  is  seen  that  the  dominant  scattering  rate  is  the  2010  process.  Unlike 
the  earlier  LO-phonon  analysis,  the  electron-electron  scattering  rates  are  not  constant  but 
rather  as  expected,  increases  with  carrier  concentration.  The  next  highest  scattering  rates 
are  for  the  2221  process,  and  the  lowest  scattering  rates  are  for  the  1110  process.  The 
total  electron-electron  scattering  rates  are  the  sum  of  the  above  processes,  taking  into 
account  the  number  of  electrons  that  change  energy  levels.  Figure  VI.A.3(a)  shows  the 
total  electron-electron  scattering  rates  for  the  2  to  1  and  the  1  to  0  transitions,  based  on 
the  processes  for  the  simplified  rate  equation  analysis.  It  is  clear  that  the  2  to  1  scattering 
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rates  are  higher  than  the  1  to  0  scattering  rates  for  all  the  initial  carrier  populations.  Thus 
the  lifetimes  are  such  that  the  electron-electron  scattering  process  by  itself  is  not 
favorable  for  the  purpose  of  keeping  a  population  inversion. 

Although  our  analysis  so  far  has  been  in  the  framework  of  the  simplified  rate 
equation  analysis  and  has  considered  only  processes  within  that  framework,  it  is  worth 
mentioning  that  in  general  a  system  with  N  number  of  energy  levels  will  have  4N  electro¬ 
electron  scattering  processes.  Some  of  these  processes  do  not  directly  by  themselves 
affect  subband  populations,  i.e.,  intrasubband  iiii  scattering  rates  where  no  net  change 
occurs  in  the  number  of  electrons  in  the  subband.  To  see  what  the  scattering  rates  are  for 
some  of  the  additional  processes,  computed  in  Figure  VI. A. 3(b)  are  the  mean  electron- 
electron  scattering  rates  for  the  additional  21 10,  2210,  and  2220  scattering  processes.  It  is 
seen  that  of  these  three,  the  2110  process  is  dominant  and  involves  a  net  change  of  losing 
one  electron  from  state  2  and  gaining  one  electron  in  state  0,  a  process  by  which  2  —*■  1 
and  1  — >  0.  The  scattering  rate  due  to  this  process  is  around  an  order  of  magnitude  less 
than  those  considered  earlier. 


Figure  VI.A.2.  Mean  electron-electron  scattering  rates  for  the  4.3  THz  step  quantum  well 
as  a  function  of  initial  state  populations  for  Teiectron  =  100  K.  The  various  scattering 
processes  for  the  2  to  1  and  the  1  to  0  transitions  for  the  simplified  rate  equation  analysis 
are  shown. 
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Figure  VI. A. 3.  Scattering  rates  for  the  4.3  THz  step  quantum  well  for  Teiectron  =  100  K.  (a) 
The  total  electron-electron  scattering  rates  for  the  2  to  1  and  the  1  to  0  transitions,  based 
on  the  processes  for  the  simplified  rate  equation  analysis,  as  a  function  of  initial  state 
populations,  (b)  The  mean  electron-electron  scattering  rates  as  a  function  of  initial  state 
populations  for  Teiectron  =  100  K,  for  some  additional  scattering  processes  for  the  2  to  1 
transition  that  are  not  within  the  framework  of  the  simplified  rate  equation  analysis. 


The  analysis  indicates  that  the  electron-phonon  scattering  rates  are  around  an 
order  of  magnitude  larger  than  the  electron-electron  scattering  rates.  The  scattering  rate 
analysis  of  the  biased  step  quantum  well  indicates  that  this  type  of  structure  is  likely 
capable  of  keeping  a  population  inversion  for  the  temperatures  and  range  of  carrier 
concentrations  considered. 

In  order  to  be  used  in  a  QC  structure,  the  step  quantum  well  must  have  an  injector 
section  to  feed  the  next  section  of  the  cascade.  One  potential  implementation  is  illustrated 
in  Figure  VI.A.4  where  one  section  is  outlined.40  The  double  barrier  injector  section  was 
designed  to  maintain  the  upper  state  lifetime,  while  allowing  the  electrons  in  the  lower 
states  to  resonantly  tunnel  to  the  upper  state  of  the  next  adjacent  section.  The  conduction 
band  profile  was  solved  by  using  a  self  consistent  solution  to  Schrodinger’s  and  Poisson’s 
equations.  The  step  well  has  been  shaped  to  spatially  separate  the  upper  state  and  lower 
states  as  much  as  possible,  in  an  effort  to  reduce  scattering  between  those  states.  The 
radiative  transition  takes  place  between  states  4  and  3,  where  E 43  =  16.6  meV  (4  THz  or  A 
~  75  //m).  The  LO-phonon  assisted  transition  takes  place  from  state  3  to  the  lower  triplet 
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states  (2,  1,  0),  since  £31  =  36.6  meV  which  is  near  hcoLo.  Scattering  between  the  triplet 
states  (2,  1,  0)  and  injection  into  the  next  adjacent  sections  upper  state  is  to  take  place. 
The  step  well  has  been  arranged  so  that  state  4  is  above  the  step  and  state  3  below  the 
highest  point,  in  an  effort  to  spatially  separate  the  wavefunctions  for  efficient  injection  to 
the  upper  state  4  and  for  reducing  the  parasitic  injection  to  state  3.  Thus,  the  step  barrier 
helps  to  reduce  unwanted  injection  into  the  middle  state,  which  is  a  problem  with  other 
square  well  LO-phonon  structures. 


Figure  VI.A.4.  Conduction  band  profde  of  a  proposed  step  well  QC  structure  found  using 
a  self  consistent  solution  to  Schrodinger’s  and  Poisson’s  equations.  The  step  well  and 
resonant  tunneling  double  barriers  of  one  section  are  outlined.  Beginning  with  the  left 
injector,  the  AlxGai  xAs  layers  compositions  are  0.143/0.035/0/0.143/0/0.143/0  and 
thicknesses  in  nm  are  4.3/20.9/13.5/1.7/9.6/2.5/7.6,  with  £43  =  16.6  meV  (4  THz)  and  £31 
=  36.6  meV.  The  applied  field  is  9.5  kV/cm.  The  9.6  nm  well  is  doped  2.9x  1 016  cnT3, 
which  corresponds  to  a  sheet  density  of  2.8X  10  cm  . 

The  mean  electron-phonon  scattering  lifetimes  for  this  structure  are 
z"3->(2,i,o)  =  0.4  psec,  for  the  upper  to  lower  states  scattering  =  1.4  psec,  and 

Z43  =  7.9  psec.  Again,  because  these  are  intrawell  transitions,  the  scattering  rates 
computed  are  relatively  fast.  The  radiative  transition  is  also  an  intrawell  transition,  and 
the  overlap  of  the  wavefunctions  is  reasonably  good,  which  gives  an  oscillator  strength  of 
about  1.34.  This  oscillator  strength  is  higher  than  the  highest  reported  LO-phonon 
designs  (~  40%  increase),  and  the  section  length  is  about  the  same,  i.e.,  about  half  that 
of  typical  bound  to  continuum  designs. 
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To  get  an  approximate  estimate  of  what  might  be  expected  for  this  type  of 
structure,  we  can  estimate  from  the  gain  equation  (equation  (III.D.27))  the  population 
inversion.  If  we  assume  a  total  loss  of  ~  40  cm-1  (which  is  reasonable  for  surface  plasmon 
waveguides)  and  assume  a  FWHM  ~  4  meV,  the  estimated  2D  population  at  threshold  is 
~  2.5x10  cm  .  This  is  within  the  carrier  densities  used  in  our  analysis.  The  current 
density  can  be  estimated  from  A/;2d  ~  (J/e)  Z4(  1  -  zV  Z43  )  and  is  calculated  to  be  J  ~  760 
A/cm-2. 

It  was  concerned  that  the  triplet  of  lower  states  spacing  was  too  close  to  the 
radiative  energy  spacing  and  could  provide  a  means  for  unwanted  absorption.  Further, 
more  energy  separation  was  wanted  between  the  lower  states  and  an  unused  state  in  the 
next  adjacent  section.  For  this  reason,  it  was  decided  to  focus  on  single  well  injector 
designs.  There  were  a  number  of  different  step  well  structures  that  were  considered,  some 
of  which  are  shown  in  Figure  VI. A. 5.  Other  structures  (not  shown),  that  resonantly 
tunneled  the  middle  and  lower  states,  step  well  designs  without  all  three  energy  levels  in 
the  same  well,  and  some  structures  that  would  also  have  required  growth  with  more  than 
two  A1  compositions  in  an  AlxGai_xAs  heterostructure,  were  considered  as  well.  Due  to 
the  difficulty  in  growing  a  sample  of  repeated  sections  with  more  than  two  A1 
compositions  in  most  MBE  chambers,  structures  requiring  only  two  A1  compositions 
were  favored. 


Figure  VI. A. 5.  Some  of  the  step  well  QC  structures  considered,  (a)  one  step,  one  well 
injector,  (b)  two  step,  one  well  injector,  and  (c)  one  step,  two  well  injector. 
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By  using  a  single  well  injector,  the  anticrossing  between  the  lower  ground  states 
can  be  kept  much  smaller  than  the  optical  transition  energy  spacing  (as  there  will  only  be 
a  doublet  of  states),  thereby  reducing  the  possibility  of  absorption.  Initially  a  single  step 
well  and  one  well  injector  (Figure  VI.A.5(a))  structure  was  analyzed.  To  avoid  a  potential 
parasitic  current  channel  (from  alignment  with  an  unwanted  energy  state  in  the  next 
adjacent  section),  the  injector  well  was  narrowed  so  that  the  only  energy  state  from  the 
well  was  that  of  the  lowest  ground  state.  This  increased  the  sub-picosecond  LO-phonon 
lifetime  of  the  middle  state  to  ~  1  psec  which  was  too  slow.  Rather  use  a  wider  injector 
well  to  keep  the  sub-picosecond  LO-phonon  scattering  lifetime  and  also  introduce  a 
second  state  from  the  well  which  could  potentially  produce  a  parasitic  current  channel,  it 
was  chosen  to  utilize  a  two  step  well  and  one  well  injector  as  shown  in  Figure  VI.A.5(b) 
(where  the  two  steps  are  of  equal  A1  composition).  This  structure,  QC  structure 
D619F10D,  will  be  discussed  in  the  next  subsection.  The  disadvantage,  however,  is  that  it 
introduces  two  additional  growth  layers  (7  growth  layers  as  compared  to  only  5  growth 
layers  for  a  single  step  well  one  well  injector  design).  Nevertheless,  even  with  a  two 
equal  A1  composition  step  well,  this  structure  is  quite  simple  having  only  two  wells. 
Albeit  the  first  being  a  step  well.  It  also  has  only  one  additional  layer  compared  to  the 
simplest  square  well  LO-phonon  QC  structure.39  Note  that  in  principle,  a  step  well  QC 
structure  does  not  need  an  injector  well  (in  which  it  would  be  a  one  well  design). 

1.  Structure  D619F10D 

The  conduction  band  profile  for  QC  structure  D619F10D  is  shown  in  Figure 
VI. A. 6.  Electrons  are  injected  into  the  upper  state  3  and  the  radiative  transition  occurs 
between  state  3  and  middle  state  2,  where  E 32  =15.1  meV  (~  3.6  THz  or  2  ~  82  /an).  The 
LO-phonon  assisted  transition  takes  place  from  the  middle  state  2  to  the  lower  doublet 
states  (1,  0),  since  £21  =  37.2  meV  which  is  near  ha>w •  This  ensures  fast  depopulation  of 
state  2  via  LO-phonon  scattering,  with  a  mean  scattering  rate  T^-qi.o)  =  0.68  psec  (0.44 
psec  band  edge).  The  anticrossing  between  states  2  and  5'  is  relatively  small  due  to  the 
thick  collector  and  injector  barriers,  and  therefore  should  not  be  a  primary  path  of 
scattering  out  of  state  2.  Though  any  electrons  that  transition  from  state  2  to  5'  would 
serve  to  help  the  population  inversion,  and  state  5'  is  near  the  resonant  LO-phonon  energy 
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spacing  to  upper  state  3'.  Injection  from  the  lower  doublet  states  (1,  0)  is  to  take  place 
into  the  next  adjacent  sections  upper  state  3’,  and  the  injector  anticrossing  Eoy  ~  1  meV. 
The  step  well  has  been  arranged  so  that  state  3  is  above  the  step  and  the  state  2  below  the 
highest  point  so  as  to  provide  spatial  separation  at  the  injector  barrier  for  improved 
injection. 


Figure  VI.A.6.  Conduction  band  profde  of  structure  D619F10D  found  using  a  self 
consistent  solution  to  Schrodinger’s  and  Poisson’s  equations.  The  step  well  and  resonant 
tunneling  well  of  one  section  are  outlined.  Beginning  with  the  left  injector,  the 
AlxGai_xAs  layers  compositions  are  0.16/0.05/0/0.05/0/0.16/0  and  thicknesses  in  nm  are 
4.8/27.9/1.8/2/6.3/4.1/6.8,  with  E^  =  15.1  meV  (3.6  THz)  and  E21  =  37.2  meV.  The 
applied  field  is  10.5  kV/cm.  The  center  2  nm  of  the  6.8  nm  well  is  doped  to  a  sheet 
density  of3.4><1010  cnT2 . 

Since  the  transitions  are  intrawell  in  nature,  the  LO-phonon  scattering  is  relatively  fast. 
The  radiative  transition  is  also  intrawell,  and  the  overlap  of  the  wavefunctions  is 
reasonably  good,  giving  an  oscillator  strength  of  -  1.03  at  9.9  kV/cm  (53.2  mV/section). 

In  order  to  determine  the  threshold  gain  required  for  lasing,  two  different 
commonly  used  waveguide  structures,  the  surface  plasmon  and  metal-metal  waveguides 
were  considered.  The  threshold  gain  was  determined  for  both  surface  plasmon  and  metal- 
metal  waveguide  resonator  configurations  at  the  operating  frequency.  The  analysis  was 
performed  assuming  an  active  region  thickness  of  -  10  jum  (which  corresponds  to  N  ~ 
185  sections)  and  a  guide  width  of  200  jum.  This  is  because  the  threshold  gain  is  worse 
(higher)  for  smaller  width  ridges.  The  threshold  gain  =  (ain+aw)IV,  where  am  is  the 
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mirror  loss  and  aw  is  the  waveguide  loss.  The  waveguide  losses  and  reflectivities  for  the 
metal-metal  waveguide  were  found  using  finite  element  method  (FEM)  solvers.108,109  The 
Drude  model  was  used  to  determine  the  material  parameters,  and  a  background  carrier 
concentration  of  2><1015  cm"1  was  assumed.  Due  to  the  difficulty  in  accurately  modeling 
the  waveguide  with  N  thinly  doped  active  region  section  layers,  the  doping  in  these  thin 
layers  was  neglected.  Thus,  the  actual  waveguide  losses  may  be  higher. 

Considering  the  surface  plasmon  waveguide  first,  the  side  contacts  spacing  was 
set  at  50  nm  to  keep  the  mode  from  coupling  to  avoid  higher  waveguide  loss.  The  top 
contact  plasma  layer  was  60  nm  thick  and  doped  at  5><1018  cnT3.  To  determine  a  suitable 
lower  plasma  layer  thickness,  the  threshold  gain  was  computed  as  a  function  of  lower 

IQ  IQ  _ T 

plasma  layer  thickness  for  doping  concentrations  of  1><10  to3x10  cm  .  The  mirror 
losses  have  been  computed  for  resonator  lengths  of  1  and  2.5  mm  with  one  facet  of  the 
waveguide  assumed  to  have  a  FIR  coating.  Figure  VI.A.7  shows  that  for  the  lowest 
threshold  gain  over  the  range  of  mirror  losses  considered  (resonator  lengths  of  1  and  2.5 

1  O  1  O  _ T 

mtn),  the  optimal  lower  plasma  layer  is  one  doped  ~  2x10  to  3X10  cm  at  a  thickness 
of  400  to  500  nm. 


Plasma  Layer  Thickness  (nm) 


Figure  VI.A.7.  Threshold  gain  for  a  200  /an  ridge,  10  /an  thick  surface  plasmon 
waveguide,  as  a  function  of  lower  plasma  layer  thickness  for  different  doping 

_ T  #  _ 1 

concentrations  (cm  )  and  mirror  losses  (cm  ). 
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The  cross-section  2D  mode  intensity  and  ID  profile  mode  intensity  are  shown  in 
Figure  VI.A.8  for  a  lower  plasma  thicknesses  of  t piasma  =  500  mn,  doped  at  3 x  1 0 1 8  cnT3. 
The  confinement  factor  is  T  =  0.28  and  the  waveguide  loss  is  aw  =  3.3  cnT1  for  this 
waveguide. 
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Figure  VI.A.8.  200  jui n  ridge,  10  /an  thick  surface  plasmon  waveguide  with  t piasma  =  500 

1  O  _ T 

mn  doped  at  3x10  cm  ,  (a)  cross-section  2D  mode  intensity,  (b)  ID  profile  mode 
intensity. 

Next,  the  threshold  gain  was  computed  for  a  metal-metal  waveguide,  where  the 
top  and  lower  contact  plasma  layers  were  60  mn  thick  (10  and  50  nm  layers  doped  at 

IQ  1  O  _ T 

5 x  1 0  and5x10  cm  ,  with  the  lower  contact  also  having  a  3.5  nm  LTG  GaAs  layer). 
For  the  metal-metal  waveguides,  we  assume  both  facets  are  uncoated.  We  find  the 
threshold  gain  to  be  7.7  and  5.7  cnT1  for  mirror  losses  corresponding  to  resonator  lengths 
of  1  and  2.5  mm  respectively.  The  attenuation  was  calculated  to  be  aw  =  4.3  cnT1  with  a 
confinement  factor  T  ~  1 .  This  near  unity  confinement  factor  is  seen  as  evident  from  the 
cross-section  2D  mode  intensity  graphed  in  Figure  VI.A.9.  Optimal  guide  parameters  for 
a  surface  plasmon  waveguide  depend  on  having  a  relatively  thick  lower  contact  layer 
compared  to  the  upper  contact  layer,  while  for  a  metal-metal  guide  both  contact  layers 
should  be  relatively  thin  to  minimize  waveguide  loss.  Because  an  etch-stop  layer  can  be 
used  for  processing  both  surface  plasmon  and  metal-metal  waveguides  from  the  same 
sample,  a  sample  that  has  an  optimal  lower  plasma  layer  thickness  for  a  surface  plasmon 
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waveguide  will  not  be  optimal  for  a  metal-metal  waveguide.  As  an  example,  the 
threshold  gain  was  computed  for  a  metal-metal  waveguide  with  a  top  contact  layer  60  nm 

1  o  _ o 

thick  doped  at  5x  10  cm  ,  and  for  lower  plasma  thicknesses  of  t piaSma  =  150  and  500  nm 
doped  at  3  x  1 0 18  cnT3  (and  found  to  change  very  little  for  guide  widths  from  100  to  200 
/an).  The  threshold  gain  was  calculated  to  be  only  8.5  cnT1  for  t piasma  =  150  nm  and 
increased  to  18.2  cnT1  for  t piasma  =  500  nm  (both  with  a  resonator  length  of  1  mm).  The 
thinner  plasma  layer  thickness  of  150  nm  can  be  seen  to  result  in  an  improvement  for  the 
metal-metal  waveguide.  Although  the  threshold  gain  is  higher  for  a  surface  plasmon 
waveguide,  the  output  coupling  factor  am/(am+au)  can  be  smaller  for  some  surface 
plasmon  waveguides  which  is  better  for  high  power  applications,  i.e.,  since  the  output 
power  =  V2r/iam/(am+aH)Nhco /\e\(I-Ith)  from  a  single  facet  of  an  uncoated  device. 


Figure  VI.A.9.  200  /an  ridge,  10  /an  thick  metal-metal  waveguide,  cross-section  2D 
mode  intensity. 


a.  MBE  Growth  and  Fabrication 

The  sample  was  grown  by  MBE  using  IQE  foundry  service,  with  N  =  185 
sections,  where  the  top  and  lower  contact  plasma  layers  were  60  nm  thick  (10  and  50  nm 
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layers  doped  at  5><10  and  5 x  1 0  cm  ,  with  the  to  be  lower  contact  also  having  a  3.5 
nm  LTG  GaAs  layer).  The  sample  was  processed  by  the  National  Research  Council 
(NRC)  (Canada)  into  metal-metal  waveguides  of  100  and  200  jui n  wide,  and  1  mm 
resonator  lengths,  with  both  facets  uncoated  (Figure  VI. A.  10).  The  fabricated  devices 
showed  an  unusually  high  impedance  (>  1  kO)  at  room  temperature,  compared  to  typical 
QCLs  of  about  10  Q.  This  likely  indicates  that  there  was  either  a  growth  problem,  such  as 
a  high  p-type  background  concentration  during  MBE  growth  (which  depleted  the  electron 
concentration  in  the  active  region),  or  non-ohmic  contacts  generated  during  the 
processing.  Initial  electroluminescence  measurements  were  made  at  T  =  10  K  (using  a 
CTI  Cryogenics  Model  22  closed  cycle  He  cryostat  and  a  Nexus  870  FT-IR  ESP  Nicolet 
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spectrometer),  and  no  light  out  was  detected.  Measurements  on  the  sample  showed  little 
current  through  the  device  near  the  injector  anticrossing  (at  -10.6  V,  -  300  mA,  which 
corresponds  to  150  A/cm  ),  as  illustrated  in  Figure  VI. A.  11  for  a  200  /an  wide  device. 
The  I-V  characteristics  were  measured  on  both  the  100  and  200  /an  wide  samples  out  to 
about  -  5  A,  without  the  device  going  into  NDR.  Because  the  contact  resistance  was  so 
high  and  no  light  was  detected,  there  was  either  a  growth  or  contact  processing  error  with 
this  sample. 


(b) 


Figure  VI.A.10.  Sample  D619F10D.  (a)  Top  view  of  the  metal-metal  waveguides, 
showing  bond  wires  to  the  100  and  200  /an  wide  by  1  mm  long  waveguides,  (b)  Cross- 
section  of  one  of  the  100  /an  wide  waveguides. 
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We  are  currently  investigating  the  cause  of  the  high  impedance  to  determine  the  origin 
which  will  help  us  refine  the  growth  and  fabrication  processes.  Tests  are  being 
performed,  using  secondary  ion  mass  spectrometry  (SIMS)  at  NRC,  to  determine  if 
correct  doping  concentrations  were  maintained  in  the  sample  growth.  No  conclusion  on 
the  feasibility  of  this  device  can  be  drawn  from  these  measurements,  as  clearly  there  was 
a  fabrication  problem. 


Figure  VI. A.  1 1.  The  measured  I-V  curve  of  sample  D619F10D,  for  a  200  //m  wide  by  1 
mm  long  metal-metal  waveguide  device. 

b.  Monte  Carlo  Simulation 

A  Monte  Carlo  simulation  of  QC  structure  D619F10D  was  performed,  to 
model  the  electron  ensemble.  As  previously  discussed,  the  scattering  between  more  than 
one  section  is  taken  into  account  because  of  the  nature  of  the  step  well  QC  structure 
(which  potentially  can  further  lead  to  an  overestimation  of  the  scattering  between  weakly 
coupled  states).  The  scattering  mechanisms  that  were  included  are  LO-phonon,  electron- 
electron,  impurity,  and  interface  roughness  scattering.  The  scattering  rates  were 
calculated  via  Fermi’s  golden  rule.  All  rates  from  these  simulations  are  net  rates  and 
include  backscattering.  State  blocking  and  screening  are  also  included  in  the  simulations. 
The  Monte  Carlo  simulations  of  the  step  well  QC  structure  were  performed  with  Tiattice  = 
25  K.  Shown  in  Figure  VI. A.  12  is  the  calculated  current  density  and  gain  (assuming  a  4.1 
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meV  FWHM)  over  the  high  gain  bias  region.  The  calculated  gain  is  higher  over  much  of 
this  region  than  the  threshold  gains  computed  earlier  for  the  waveguide  resonator 
configurations  analyzed.  However,  just  beyond  50  mV/section  a  dip  in  the  gain  curve  is 
observed.  It  was  felt  that  a  more  unifonn  gain  could  be  realized,  by  modifying  the 
structure  slightly  to  improve  the  injection  efficiency. 


Figure  VI. A.  12.  Results  of  the  Monte  Carlo  simulation  of  structure  D619F10D,  with  and 
without  including  interface  roughness  scattering,  (a)  Current  density  as  a  function  of 
applied  electric  field  bias,  (b)  Gain  as  a  function  of  bias. 

2.  Structure  D619F10E 

The  dip  in  the  gain  curve  seen  from  the  Monte  Carlo  simulation  of  QC  structure 
D619F10D,  was  likely  a  result  of  inefficient  scattering  into  the  upper  state  at  those  biases. 
In  order  to  correct  for  this,  a  second  structure  was  designed  and  analyzed,  structure 
D619F10E  where  the  injector  well  was  thinned  slightly  and  the  collector  well  thickness 
increased  slightly.  With  these  changes,  it  was  anticipated  a  more  uniform  gain  over  these 
biases  could  be  achieved.  The  conduction  band  profile  for  the  step  well  QC  structure  is 
shown  in  Figure  VI.A.13.113  Electrons  are  injected  into  the  upper  state  3  and  the  radiative 
transition  occurs  between  state  3  and  middle  state  2,  where  £32  =15.2  meV  (~  3.7  THz  or 
2,-82  //m).  The  LO-phonon  assisted  transition  takes  place  from  the  middle  state  2  to  the 
lower  doublet  states  (1,  0),  since  £21  =  37.9  meV  which  is  near  hcoio  (~  36  meV  in 
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GaAs).  This  ensures  fast  depopulation  of  state  2  via  LO-phonon  scattering,  with  a 
scattering  rate  lifetime  Z2^o,o)  ~  0.5  psec  (band  edge)  for  this  structure.  Again,  the 
anticrossing  between  states  2  and  5’  is  relatively  small  (less  than  ~  0.5  meV)  due  to  the 
thick  collector  and  injector  barriers,  and  therefore  should  not  be  a  primary  path  of 
scattering  out  of  state  2.  Though  any  electrons  that  transition  from  state  2  to  5'  would 
serve  to  help  the  population  inversion,  and  state  5'  is  near  the  resonant  LO-phonon  energy 
spacing  to  upper  state  3'.  Injection  from  the  lower  doublet  states  (1,  0)  is  to  take  place 
into  the  next  adjacent  sections  upper  state  3',  and  the  injector  anticrossing  has  been 
slightly  increased  to  E$y  ~1.2  meV.  The  step  well  has  been  arranged  so  that  state  3  is 
above  the  step  and  the  state  2  below  the  highest  point  so  as  to  provide  spatial  separation 
at  the  injector  barrier  for  improved  injection. 


Figure  VI. A.  13.  Conduction  band  profile  of  structure  D619F10E  found  using  a  self 
consistent  solution  to  Schrodinger’s  and  Poisson’s  equations,  shown  at  an  applied  electric 
field  of  10.4  kV/cm  (56.3  mV/section).  The  step  well  and  resonant  tunneling  well  of  one 
section  are  outlined.  Beginning  with  the  left  injector,  the  AlxGai  xAs  layers  compositions 
are  0.16/0.05/0/0.05/0/0.16/0  and  thicknesses  in  mn  are  4.6/27.9/1.8/2/6.3/4.7/6.8,  with 
E32  =  15.2  meV  (~  3.7  THz)  and  E21  =  37.9  meV.  The  center  2  mn  of  the  6.8  mn  well  is 
doped  to  a  sheet  density  of  3.4><1010  cm  2. 

Since  the  transitions  are  intrawell  in  nature,  the  LO-phonon  scattering  is  relatively  fast. 
The  radiative  transition  is  also  intrawell,  and  the  overlap  of  the  wavefunctions  is 
reasonably  good,  giving  an  oscillator  strength  of  ~  0.94  at  9.9  kV/cm  (53.6  mV/section). 

The  Monte  Carlo  method  was  used  to  model  the  electron  ensemble  of  QC 
structure  D619F10E,  in  the  same  manner  as  for  the  previous  structure  analyzed.113  The 
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Monte  Carlo  simulations  of  the  step  well  QC  structure  in  Figure  VI.A.13  were  performed 
with  T lattice  =  25  K.  Table  VI. A.I  lists  the  subband  populations  and  electron  temperatures 
found  at  an  applied  electric  field  of  53.6  mV/section. 


n 

«2D  (xlO9  cm  2) 

Telectron  (K) 

0 

16.1 

131 

1 

10.1 

132 

2 

1.53 

155 

3 

6.22 

127 

Table  VI.  A.I.  Populations  of  the  subbands  ( En )  and  electron  temperatures  found  from  the 
Monte  Carlo  simulation  of  structure  D619F10E,  at  an  applied  electric  field  of  53.6 
mV/section,  with  Tiattice  =  25  K. 

The  current  density,  gain,  and  population  densities  and  the  electron  temperatures  for  the 
upper  (£3)  and  middle  (E2)  states  are  shown  in  Figure  VI.A.14  over  a  range  of  applied 
electric  field  biases.  The  results  are  shown  with  and  without  including  interface 
roughness  scattering.  Others  have  shown  reasonable  agreement  with  experimental 
findings  by  not  including  interface  roughness  scattering  in  their  simulations,114  and  that 
agrees  with  our  simulations  of  other  square  well  structures.  In  these  simulations  interface 
roughness  scattering  was  found  to  increase  somewhat  the  current  density  and  electron 
temperatures,  and  also  slightly  decreases  the  medium  gain.  The  current  density  as  a 
function  of  applied  bias  follows  about  what  is  expected,  with  the  exception  of  the  points 
between  ~  45  to  50  mV/section.  This  is  where  some  scattering  to  parasitic  current 
channels  takes  place.  These  high  values  are  likely  due  to  an  overestimation  of  the 
scattering  to  the  weakly  coupled  states.  A  small  bump  in  the  curve  at  ~  38.4  mV/section 
is  observed,  where  this  is  near  the  anticrossing  between  states  0  and  2’.  The  step  serves  to 
also  minimize  scattering  between  states  0  and  2'  at  this  intermediate  bias. 


93 


Figure  VI.A.14.  Results  of  the  Monte  Carlo  simulation  of  structure  D619F10E,  with  and 
without  including  interface  roughness  scattering,  (a)  Current  density  as  a  function  of 
applied  electric  field  bias,  (b)  Gain  as  a  function  of  bias,  (c)  The  population  densities  and 
(d)  electron  temperatures  for  the  upper  (C3)  and  middle  (E2)  subbands  as  a  function  of 
bias. 


A  peak  gain  of  86.9  cnT1  is  found  (assuming  a  FWHM  =4.1  meV)  and  the  gain  is 
~  70  cm1  or  higher  across  much  of  the  high  gain  bias  region.  This  is  due  to  the  vertical 
transitions  of  the  optical  and  LO-phonon  transitions,  as  well  as  the  high  injection 
efficiency  from  the  step.  The  vertical  optical  transition  allows  for  a  relatively  high 
oscillator  strength  and  fast  LO-phonon  scattering  helps  keep  the  middle  state  2  depleted. 
Due  to  the  step,  injection  loss  to  the  middle  state  2'  is  relatively  small.  The  injection 
efficiency  is  found  to  be  ~  89%,  defined  here  to  be  from  the  lower  doublet  states  (0,  1) 
to  the  upper  state  3'.  The  lifetime  of  the  middle  state  due  to  LO-phonon  scattering  is 
found  to  be  t2->( \ ,o)  =  0.88  psec,  with  the  ratio  of  from  LO-phonon 

scattering  varying  from  about  2.1  to  greater  than  3  across  biases  in  the  high  gain  region. 
The  calculated  gain  is  seen  to  be  higher  than  the  calculated  threshold  gains  for  both  the 
surface  plasmon  and  metal-metal  waveguides  analyzed  earlier. 
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B.  SQUARE  WELL  LO-PHONON  QC  STRUCTURE 

As  previously  discussed,  only  square  quantum  well  structures  have  been  used  in 
QCL  designs.  The  highest  power  QCL,  reported  in  the  given  reference,10  utilized  a  four 
square  well  LO-phonon  structure  and  operated  at  ~  4.4  THz.  This  device  represents  a 
relatively  optimized  design  in  terms  of  the  structure,  wavelength  of  operation,  and 
waveguide,  to  achieve  a  peak  power  of  ~  250  mW  in  a  ~  200  /an  wide  by  ~  1.2  mm  long 
surface  plasmon  waveguide  configuration.  There  were  N  =  183  sections  grown,  and  the 
lower  plasma  layer  was  400  mn  doped  at  3><10  cm  .  As  with  all  conventional  square 
well  LO-phonon  designs,  it  contained  the  usual  first  three  square  wells,  and  a  fourth 
injector  well  similar  to  the  first  LO-phonon  THz  QCL.  The  third  well  is  used  to 
resonantly  tunnel  the  lower  lasing  state  and  also  arrange  for  a  lower  ground  state. 
Because  of  this,  as  with  other  LO-phonon  square  well  designs,  there  is  a  doublet  of  states 
at  the  lower  lasing  state.  This  structure  will  next  be  described  in  detail,  and  we  will  model 
it  using  our  Monte  Carlo  code  to  compare  the  simulation  results  to  the  experimentally 
measured  results. 

1.  4.4  THz  Structure 

The  conduction  band  profile  for  the  4.4  THz  QC  structure  is  shown  in  Figure 
VI.B.l.  Electrons  are  injected  into  the  upper  state  4  and  the  radiative  transition  occurs 
between  state  4  and  state  3.  The  design  was  intended  to  emit  near  4.7  THz,  but  was 
experimentally  found  closer  to  ~  4.4  THz  (k  ~  68  //in).  The  LO-phonon  assisted 
transition  primarily  takes  place  from  state  2  (as  well  as  from  state  3)  to  the  lower  doublet 
states  (1,0),  since  £21  is  near  ha>w  (~  36  meV  in  GaAs).  States  3  and  2  are  intended  to  be 
coupled  for  fast  scattering  between  those  states.  This  ensures  fast  depopulation  of  the 
lower  lasing  state  via  LO-phonon  scattering.  Injection  from  the  lower  doublet  states  (1,0) 
is  to  take  place  into  the  next  adjacent  sections  upper  state  4'.  This  shorter  wavelength  of 
operation  and  hence  larger  radiative  energy  separation  of  the  upper  and  lower  lasing 
states,  along  with  good  coupling  at  the  injector,  aids  in  improved  injection. 
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Figure  VI.B.l.  Conduction  band  profile  of  the  4.4  THz  structure,  shown  at  an  applied 
electric  field  of  1 1.5  kV/cm  (62.4  mV/section),  with  one  section  outlined.  Beginning  with 
the  left  injector,  the  AlxGai-  xAs  layers  compositions  are  0.15/0/0.15/0/0.15/0/0.15/0  and 
thicknesses  in  nm  are  4.8/8.2/1.7/6.8/4.0/16.4/3.4/9.0.  The  16.4  nm  well  is  doped  to  a 
sheet  density  of  3.1X1010  cnT2. 

The  Monte  Carlo  simulations  of  the  4.4  THz  QC  structure  were  performed  with 
T lattice  =  25  K.  The  current  density  and  gain  are  shown  in  Figure  VI. B. 2  over  the  high  gain 
bias  region  of  the  device.  The  results  are  shown  with  and  without  including  interface 
roughness  scattering,  for  roughness  of  all  (same  as  previously  used  for  the  step  well 
structures  analysis)  and  a/4.  In  these  simulations,  reasonable  agreement  is  shown  with  the 
interpolated  measured  values10  for  the  current  density  at  most  biases,  with  the  largest 
difference  being  an  overestimation  occurring  ~  68  mV/section  which  is  just  beyond  the 
injector  anticrossing  and  near  the  anticrossing.  This  overestimation  is  likely  due  to 
not  including  coherent  effects  of  dephasing.  As  with  the  previous  step  well  structure 
analyzed,  the  effects  of  interface  roughness  scattering  were  found  to  increase  somewhat 
the  current  density,  and  also  slightly  decrease  the  medium  gain.  It  is  seen  that  assuming  a 
roughness  of  all  overestimates  the  current  density,  and  that  using  a  roughness  of  a/4  is 
almost  negligible  compared  to  simulations  without  including  interface  roughness 
scattering.  Although  one  monolayer  (ML)  is  sometimes  listed  as  a/2,  this  really  depends 
on  how  it  is  defined.  A  GaAs  unit  cell  consists  of  two  Ga  layers  and  two  A1  layers. 
During  growth,  it  is  possible  to  monitor  the  Ga  and  As  layers  separately.  If  each  of  these 
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layers  is  defined  as  one  monolayer,  then  there  are  four  monolayers  and  1  ML  =  a/4  by 
definition.  On  the  other  hand,  if  one  defines  one  Ga  layer  plus  one  As  layer  as  one 
monolayer,  then  there  are  two  monolayers  and  1  ML  =  a/2  by  definition.  The  gain  is  also 
calculated,  and  the  device  yields  high  gain. 


Figure  VI.B.2.  Results  of  the  Monte  Carlo  simulation  of  the  4.4  THz  structure,  with  and 
without  including  interface  roughness  scattering,  (a)  Current  density  as  a  function  of 
applied  electric  field  bias,  (b)  Gain  as  a  function  of  bias. 

This  analysis  indicates  that  interface  roughness  scattering  effects  are  only 
appreciable  for  roughness  greater  than  ~  1  ML,  and  that  overestimation  in  the  current 
density  may  occur  assuming  a  roughness  of  a/2.  The  results  presented  earlier  for  the  step 
well  structure,  which  included  interface  roughness  scattering  assuming  a  roughness  of 
a/2,  are  likely  to  be  overestimating  the  current  density  of  the  structure.  While  interface 
roughness  can  vary  from  sample  to  sample,  this  analysis  shows  the  importance  of 
maintaining  correct  tolerances  during  growth.  Because  our  analysis  shows  reasonable 
agreement  with  the  experimentally  measured  current  densities,  it  indicates  that  our  Monte 
Carlo  analysis  of  the  step  well  structures  should  have  similar  accuracy. 
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VII.  CONCLUSIONS  AND  DISCUSSION 


This  research  focused  on  modeling  electron  transport  properties  of  THz  LO- 
phonon  QC  structures,  including  the  design  and  analysis  of  new  step  well  QC  structures. 
From  the  Monte  Carlo  simulations,  it  was  found  that  step  well  structures  are  capable  of 
high  gains  and  oscillator  strengths,  high  injection  efficiencies,  with  comparable 
characteristics  to  other  square  well  designs,  but  do  have  increased  scattering  from  the 
upper  state  to  the  lower  sates.  The  peak  gain  of  ~  87  cm  1  found  in  these  simulations  of  a 
step  well  QC  structure,  is  higher  than  the  peak  gain  of  ~  73  cm  1  reported  for  a  simulated 
LO-phonon  square  well  design,  which  also  showed  similar  overestimation  of  the  current 
density  from  their  Monte  Carlo  simulations.114  The  current  density  near  the  injector 
anticrossing  was  higher  in  these  step  well  QC  structure  results  and  some  of  the  electron 
temperatures  were  slightly  higher  as  well.  This  overestimation  of  the  current  density  was 
likely  due  to  not  taking  into  account  coherence  effects. 

Perhaps  one  of  the  most  important  results  from  these  simulations  on  the  step  well 
QC  structures,  showed  high  injection  efficiencies  are  possible  using  a  step  well  injector. 
Though  the  injection  efficiency  was  much  higher  than  typical  for  square  well  LO-phonon 
designs  with  similar  oscillator  strengths  (~  1 .6  times  higher),  the  current  density  was  still 
high  due  to  scattering  from  the  upper  to  lower  states.  It  may  be  possible  to  utilize  a  step 
well  injector  to  improve  injection  efficiency,  even  if  the  step  is  not  used  to  arrange  all 
three  energy  states  within  the  same  well.  Utilizing  a  step  well  injector  with  resonant 
tunneling,  could  yield  the  benefits  of  both  good  injection  efficiency  and  improved  upper 
to  lower  state  lifetime.  Our  analysis  is  the  first  to  suggest  using  more  than  one  potential 
height,  using  two  or  more  A1  compositions,  in  a  QCL  active  region  structure. 

The  step  well  structure  that  was  grown  and  processed,  failed  experimentally. 
Because  of  the  large  contact  resistance  and  the  low  current  through  the  device  at  the 
injector  anticrossing,  it  was  determined  this  was  due  to  a  growth  or  contact  fabrication 
problem.  Therefore,  it  was  not  possible  to  conclude  the  practical  feasibility  of  such 
structures.  Testing  is  being  conducted  to  detennine  the  source  of  the  fabrication  problem. 
The  step  well  design  remains  a  promising  approach  that  differs  from  the  two  successful 
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approaches  that  have  previously  been  used  (miniband  and  square  well  LO-phonon). 
Because  the  middle  state  (upper  phonon  or  lower  lasing  state)  is  a  single  energy  state, 
contrasting  to  previous  LO-phonon  based  QCL  designs  that  have  doublet  states,  for  long 
wavelength  lasers  this  could  be  important  because  it  eliminates  the  possibility  for 
unwanted  THz  absorption  that  could  otherwise  occur  between  those  doublet  states. 

The  Monte  Carlo  simulations  perfonned  on  a  high  power  conventional  square 
well  LO-phonon  design,  showed  reasonable  agreement  with  experimentally  measured 
results.  This  indicated  that  our  simulations  of  the  step  well  structures,  was  likely  of 
similar  accuracy.  It  was  also  found  that  interface  roughness  scattering  becomes 
significant  for  roughness  greater  than  approximately  one  monolayer,  and  including 
interface  roughness  scattering  in  simulations  can  lead  to  an  overestimation  of  the  current 
density  in  the  device,  depending  on  the  roughness  height  chosen.  This  illustrates  the 
importance  of  maintaining  tight  growth  tolerances  for  QC  structures. 

Although  a  high  oscillator  strength  can  lead  to  high  gain,  provided  the  upper  state 
lifetime  can  be  maintained,  optimum  temperature  devices  may  be  realized  using  lower 
oscillator  strength  structures.  It  should  be  pointed  out  that  the  best  temperature 
performing  devices  are  due  to  their  higher  frequency  of  operation,  and  none  across  the 
THz  frequency  range  perform  much  better  than  T  ~  EradiativelkB  =  hco/kBn'l~  A  longer 
wavelength  device  might  operate  at  a  lower  temperature,  but  with  respect  to  EmdiativJkB  be 
the  same  or  better  than  a  shorter  wavelength  device  of  higher  operating  temperature. 
Even  in  step  well  structures,  it  is  possible  to  make  the  radiative  transition  more  diagonal 
by  utilizing  an  additional  well  prior  to  the  step,  and  also  by  not  arranging  all  three  energy 
levels  within  the  same  well. 

The  overestimation  between  weakly  coupled  states  that  occurs  in  the  Monte  Carlo 
simulations,  can  be  dealt  with  either  by  including  dephasing  in  a  density  matrix  Monte 
Carlo  simulation  or  by  using  nonequilibrium  Green’s  functions.  The  nonequilibrium 
Green’s  function  technique  is  more  intuitive  from  a  quantum  view  point,  and  has  several 
advantages.  It  can  be  used  to  include  dephasing,  many-body  effects,  and  allow  for 
spectral  estimation.115'116  This  technique  is  useful  in  many  other  areas  of  condensed 
matter  physics  as  well.  We  are  continuing  to  research  these  techniques  and  structures 
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discussed,  and  are  also  looking  to  include  the  effects  of  other  parameters  such  as 
temperature  in  our  analysis. 
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